53

Finite Element: volumetric integrals, internal polynomial order Classical finite element methods assume continuous or weakly continuous approximation spaces and ask for volumetric integrals of the weak form to be satisfied. The order of accuracy is increased by raising the approximation order within elements. The methods are not exactly conservative, thus ...


19

You should definitely check out Julia. Julia is a programming language which is similar to Python or MATLAB but utilizes a strong type-inference algorithm + JIT in order to optimize code. If types can be fully inferred in a function (which it usually can), then the code compiles in a fully static manner that matches C or Fortran performance. Translating ...


12

Non-degenerate diffusion will prevent discontinuities in the true solution, but those may be poorly resolved. I will assume that your question is about how to evaluate diffusive fluxes of the form $f(u,\nabla u,p) = \kappa(u,\nabla u,p)\nabla u$ where $\kappa$ is a $d\times d$ SPD tensor and $p$ represents some spatially distributed parameters that are cell-...


12

No, positive definiteness (and symmetry) are only precondition to using the Conjugate Gradient method. But there are plenty of other iterative methods such as MinRes and GMRES that can be used for indefinite and non-symmetric matrices.


11

What is it you want to achieve? If you want to develop methods/algorithms you might prefer a language that is flexible, and that you are familiar with. As you stated in your question, the Fortran code of your professor was hard to grasp for you, so you re-implemented it in MATLAB. This is very natural way of doing method development: first you implement ...


11

Start simple. Learn Python. I have been paid to write programms for over forty years and I have used all the languages mentioned in other answers (except Julia - I had never heard of it before now). Each language has its strengths, and most have their weaknesses. Like human languages, code is a way of expressing and framing ideas and when you know how to ...


10

I was pondering this a few days ago (also in Python). Personally I don't think that object oriented programming is always a good fit for numerical programming. You can get distracted with designing the classes rather than just solving the equations. I prefer to stay with simple functions, and with numpy you can have your equations vectorised so the number of ...


9

This is rather a general remark on FVM than an answer to the concrete questions. And the message is that there shouldn't be the need for such an adhoc discretization of the boundary conditions. Unlike in FE- or FD-methods, where the starting point is a discrete ansatz for the solution, the FVM approach leaves the solution untouched (at first) but averages ...


9

The structure you show is a common choice and equivalent to storing cell-face adjacencies in a CSR matrix format, with the boundary ghost cells in a special place. However, note that FV methods can also be formulated to consist entirely or almost-entirely of face traversal where each face is only visited once (reconstruct to face centroid/quadrature point ...


9

Simple answer: in modern python every data type is a class, so formally there is no difference between the two solutions you proposed. (Please remember to use new-style classes: classic classes are obsolete! See http://docs.python.org/2/reference/datamodel.html#new-style-and-classic-classes) Now the question should be: how do I organize an efficient data ...


9

Just as an aside, your github documentation is fantastic. This is just a guess from DG methods, which can have similar issues if numerical fluxes aren't chosen carefully (I figure FV methods are a subset of DG methods). If you're using interpolation from cell centers to define your fluxes, then this should be equivalent to using the average as a numerical ...


9

The finite element method is actually quite widely used in fluid flow problems, for example for the Stokes and Navier-Stokes equations. The delineation between the methods is more along the following line: The finite element method is very well suited for second order (in space) differential equations. That has something to do with the fact that trial and ...


8

The conceptual differences of FEM and FVM are as subtle as the differences between a tree and a pine. If you compare a certain FEM scheme to the FVM discretization applied to a particular problem then you can speak of fundamental differences that become evident in different implementation approaches and different approximation properties (as @Jed Brown has ...


8

I don't really have much more to say than I already did on the pages you linked to, but for me the primary arguments go like this: In many problems, one needs to adapt the mesh between time steps. The conceptual framework for doing this is the Rothe method where one can choose the spatial discretization independently at each time step whereas the method of ...


8

The difference between finite volumes and finite differences is really more about the form of the equations solved. In typical FV methods, the conservative form is discretised in terms of integrals and fluxes, whereas FD methods generally approximate derivatives in the non-conservative form directly. Maintaining conservation and preserving physical ...


7

We have the following problem: $$\frac{\partial u}{\partial t}+v\color{red}{\frac{\partial u}{\partial x}}-\nu\color{blue}{\frac{\partial^2u}{\partial x^2}}=0 \tag{*}$$ The function $u$ may represent for example the concentration that propagates at velocity $v>0$ and disperses in a medium with viscosity $\nu>0$. Since only we are discussing how terms ...


7

I highly recommend to anyone, regardless of background, learning both low-level, "fast" language (C, Rust, C++, Go) high-level scripting language (Python, MATLAB, Mathematica, R, bash) As a general statement, lower-level languages far more precisely control hardware and are very efficient, while higher-level languages are easier to write/describe logic in. ...


6

You need to solve a Riemann problem, perhaps approximately. For a linear system of equations, the solution to the Riemann problem is just upwinding applied to the characteristics. An "exact" Riemann solver for nonlinear problems resolves the full wave structure (consisting of shocks, rarefactions, and possibly linearly degenerate contact discontinuities). An ...


6

Of course you can discretize the two equations of your system with two different methods. The challenge will simply be when the solution of one equation enters that of another. At that point, you will have to decide what the finite volume solution should be at a finite difference point should be, or how to integrate the finite difference solution (which is ...


6

The complexity of simplicial tessellations in $\mathbb{R}^{d}$, for $d>2$ is, unfortunately, not linear. In general, a simplicial tessellation $\mathcal{T}(X)$ for $n$ vertices $X\subset\mathbb{R}^{d}$ can have $\Theta(n^{\lceil d/2\rceil})$ d-simplexes in the worst case. Similar results exist for the lower dimensional faces of $\mathcal{T}$. Based on ...


6

The restrictive CFL of DG schemes typically comes from the combination of high order accuracy and a compact stencil (see this reference for example). The CFL depends on bounding the variational form in terms of the $L^2$ norm of the solution, which depends on derivative and traces of polynomials. Bounds for each of these quantities (using Bernstein or ...


6

The issue with finite elements in the current context is not that you have a first order differential equation, but with the kind of first order equation you have. In general, finite element methods can be made to perform quite well for first order problems such as advection equations. But that's beside the point here. To illustrate the issue, think of the ...


6

Of course, the Fourier transform is a linear operator. So, you have the kinetic energy defined as: $E(\mathbf{r}) = \frac{1}{2} \mathbf{u}(\mathbf{r}) \cdot \mathbf{u}(\mathbf{r})$. The Fourier transform of $\mathbf{u}$ and $E$ are: $$\tilde{E}(\mathbf{k}) = \int_{\Omega} E(\mathbf{r}) e^{-i \mathbf{k} \cdot \mathbf{r}} d^{3} \mathbf{r}$$ $$\tilde{\mathbf{...


6

You have some great answers already. I think there is no single answer to your question. What language(s) you choose to learn depends on what you intend to do. When I was a graduate student I too learned some Matlab, Mathematica, Maple, SAS, Stata, because my university had licenses and because in my area of research that's what others used, so I could ...


5

I suggest the following simple algorithm: Depending on whether you want symmetry of the mesh, either start at a point $x_0$ at the left end of your interval, or in the center. In the latter case, just mirror the mesh you obtain. Having performed the algorithm to $x_k$, compute $x_{k+1}$ as $x_k+h_k$, where \begin{gather} h_k = \min \left\{ h_{\max}, \frac{\...


5

The basic difference is simply the meaning to be attached to the results. FDM predicts point values of any aspect of the solution. Interpolation between these values is often left to the imagination of the user. FVM predicts averages of conserved variables within specific control volumes. Hence it predicts the integrated conserved variables and can be shown ...


5

Others have said it all, but I just wanted to add a simple, yet sometimes subtle, point. Your upwind discretization remains conservative as long as you use a consistent interpolation of $a(x)$ on the cell boundaries. What I mean by consistent is that the only condition that the interpolation needs to satisfy is $$ a_{i+1/2}^+ = a_{i+1/2}^- $$ In other ...


5

The finite element method has similar problems to FD with regards to stability in solving the advection-diffusion equation, i.e. the same restrictions on Peclet number apply. One remedy, also similar to that used in FD, is to use a formulation that includes "upwinding". A nice set of lecture notes that discusses how upwinding is added to FE formulations is ...


5

There are whole books written on this, but you should investigate (search the web, really) upwind diffusion and Streamline Upwind Petrov Galerkin finite element methods first. There many more methods that extend these ideas or approach the problem from a different direction.


5

The sources you are looking at are all looking at hyperbolic problems. The issues are different for elliptic problems and "limiters" are generally not the preferred tool. I outlined some of the methods and tradeoffs in this answer. As for time integration, $L$-stability is the important property to prevent bad overshoots for parabolic systems. ...


Only top voted, non community-wiki answers of a minimum length are eligible