We’re rewarding the question askers & reputations are being recalculated! Read more.
49

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 ...


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-...


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

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

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

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 ...


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 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

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 ...


6

There are no such operators so that $V2C(C2V)=I$ and $C2V(V2C)=I$ simultaneously. This is already easy to see if you only have a 1d situation with 2 vertices and 1 cell. In that case, $C2V$ is a $2 \times 1$ matrix, and $V2C$ is a $1 \times 2$ matrix. It is easy to verify that you can't find entries for these two matrices that satisfy the criteria you ask. ...


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

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 ...


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

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 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. ...


5

As Jed says, limiters are not usually an efficient approach for parabolic/elliptic problems. WENO is much more expensive than simple piecewise-polynomial interpolation, so I would first try vanilla interpolation and see if you actually have oscillations. WENO is really designed for situations in which the solution is discontinuous; yours is not. In case ...


5

Some thoughts from someone who has worked a fair amount in compiled languages, and has done a tiny bit of FVM: Typically, if you have experience programming in C, you sketch out a high-level description (pseudocode) of what you would like to do. Then you look for libraries that might implement the data structures and capabilities you need for your high-...


5

If you wish to know only the change in density you are interested in the dilatation or volumetric strain of the stress tensor associated with the transformation. If you also wish to know about the anisotropy of the stretching you should look at computing the right Cauchy-Green deformation tensor, $\Delta$, and it's associated eigenvalues and eigenvectors. ...


5

The concept of $p$-nonconforming meshes can also apply to continuous FEM, not just DG. For continuous FEM, continuity is enforced by enforcing a single set of degrees of freedom on a face shared between two elements. If those elements have varying polynomial degrees, the trace space on the face must be made the same. This may be done by restricting the ...


5

"Spectral methods" usually means methods which make use of global basis functions. Fourier spectral methods use sine and cosine, and are used when you have periodic boundary conditions. Chebyshev methods use Chebyshev polynomials and are useful in non-periodic cases. These two methods are used in DNS, see e.g., hit3d which uses fourier and periodic bc, and ...


5

Taking the average will certainly work, but is recommended to take the harmonic mean to account for the changing conductivity: $k = \Big( \frac{1-f}{k_i} + \frac{f}{k_{i+1}} \Big)^{-1}$ If you choose your faces directly midway between your two cells, f = 0.5: $ k = \frac{2k_ik_{i+1}}{k_i + k_{i+1}} $ It can be derived by looking at the physical flux, ...


5

You can't. Nonlinear systems of equations are in general not solvable exactly. What you need to do is to use a method to solve nonlinear systems, of which there are of course quite a lot: A simple approach would be to use $D(U)\approx D(U^{n-1})$, where $D^{n-1}$ is the solution of the previous time step. A possibly smarter approach would be to use $D(U)\...


5

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{...


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