24

Nitsche's method is related to discontinuous Galerkin methods (indeed, as Wolfgang points out, it is a precursor to these methods), and can be derived in a similar fashion. Let's consider the simplest problem, Poisson's equation: $$ \left\{\begin{aligned} -\Delta u &= f \qquad\text{on }\Omega,\\ u &= g \qquad\text{on }\partial\Omega. \end{aligned}\...


18

Short answer: No, you don't have to do integration for certain FEMs. But in your case, you have to do that. Long answer: Let's say $u_h$ is the finite element solution. If you choose piecewise linear polynomial as your basis, then taking $\Delta$ on it will give you an order 1 distribution (think taking derivative on a Heaviside step function), and the ...


17

Nothing stops you from doing that technically, but when you integrate by parts you get more flexibility with the solution space in that they need not have $H^2$ regularity (required for the non I.B.P formulation). The linear elements you suggest generally have enforced continuity between elements, and so could not be in $H^2$. The I.B.P formulation ...


15

Excellent answers already on this page, but there is still a (small) missing point. The OP asked: Now, let's say that I have a PDE with higher order derivatives, does that mean that there are many possible variational forms, depending on how I use Green's formula? And they all lead to (different) FEM approximations? Integrating by parts (in the ...


11

Yes, this is the standard Aubin-Nitsche (or duality) trick. The idea is to use the fact that $L^2$ is its own dual space to write the $L^2$-norm as an operator norm $$\|u\|_{L^2} = \sup_{\phi\in L^2\setminus\{0\}} \frac{(u,\phi)}{\|\phi\|_{L^2}}.$$ We thus have to estimate $(u-u_h,\phi)$ for arbitrary $\phi\in L^2$. To do that, we "lift" $u-u_h$ to $H^1_0$ ...


10

From a numerical perspective, it's perhaps easiest to discuss the discretizations directly. For the Poisson equation with homogeneous Dirichlet boundary conditions, there is a unique solution for any right-hand side. Once discretized, the equation can be written in the form $Ax = b$, where $A$ is the standard discretization of the 3D Laplacian operator with ...


9

The biharmonic equation is the Euler-Lagrange equation of the Laplacian energy $\frac{1}{2} \langle \Delta u,\Delta u \rangle$. A systematic approach to discretize higher order problems is to convert the unconstrained problem to a constrained problem: Minimize $\frac{1}{2} \langle v,v\rangle$ s.t. $\Delta u=v$; that is, \begin{equation} \frac{1}{2} \langle v,...


6

Preliminary remarks: The natural jump condition is on the normal derivative, not on the full gradient. The fact that the solution is continuous across the interface is part of the formulation of the problem. This is not only due to the finite element formulation. Formulation The standard formulation simply enforces it with... no extra term. Start by ...


6

I know its been 2 years.. but here is an example from physics to hammer it home. Consider you are solving the Poisson equation: $$\nabla^2\phi =f$$for electrostatic potentials(ϕ) with the right hand side: $$ f=(-ρ/εo) $$ And, ρ(the charge density) is specified at all the grid points. Statement: For Periodic solution, the ∫ρdv=0 condition (described by Ben ...


6

You are solving Poisson's equation: $$\nabla^2 \Phi = 4\pi G \rho.$$ Notice that if $\Phi$ is a solution, then so is $\Phi+C$ for any constant $C$. Furthermore, $C$ will have no effect on your simulation since the forces depend only on derivatives of $\Phi$. The divide-by-zero issue is just a manifestation of this degeneracy. In the Fourier-...


6

The formulation of this problem is tricky. Here is what you have in your original post: Find $h \in L^2(\partial D)$ such that for any $w \in H^{\frac 1 2}(\partial D)$, $$ \int_{\partial D} hw \,dl = \int_D \nabla w \nabla u dS. $$ This already can't be quite right because you have $w$ on the right in a domain integral, so it cannot be a function $w \in ...


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

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

For the particular equation you are solving (called the minimal surface equation), the functional you are trying to minimize is $$ J(u) = \int_\Omega \sqrt{1 + |\nabla u|^2} \; dx. $$ You can find a derivation of the equations, as well as a discussion of solution approaches in lectures 31.5 and following here: http://www.math.tamu.edu/~bangerth/videos....


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


4

Your derivation is correct up to minor imprecisions: \begin{multline} \int_\Omega \bigg((\lambda +2\mu) \nabla^2 \mathbf u \cdot \mathbf v +(\mu+\lambda)\nabla \times (\nabla \times\mathbf u)\cdot \mathbf v + \mathbf f \cdot \mathbf v \bigg) d\Omega = \\ = -\int_\Omega \bigg((\lambda +2\mu) \nabla \mathbf u :\nabla \mathbf v +(\mu+\lambda)(\nabla \times \...


4

Look for the book Vibration of plates by Arthur Leissa. It has explicit solutions for square and circular plates. Including tables with approximated eigenvalue for different boundary conditions.


4

See NCLab, which has a finite element analysis module using the hp-FEM library Hermes. For example, see the elasticity module.


4

I decided to expand my earlier comment into an answer. I'd suggest using the Morley element that uses $P_2$ basis in each element and the degrees of freedom are values at vertices normal derivatives at edge midpoints It is probably the easiest element for biharmonic problem implementation wise. Here is an example using sp.fem library (version 8019aa7): ...


4

This would call for adaptive finite element methods, and if you are only interested in specific regions, for the use of goal oriented error estimators to drive the adaptive mesh refinement. The problem you want to solve is a pretty standard one and you can find my own contribution to this in the step-6 tutorial program of the deal.II library for the general ...


4

After a many days of discussions, your problems are 1) you choose too simple problem for which you reach the machine zero and therefore cannot observe normal properties 2) perhaps you do not use the boundary values correctly. With the scheme you work the values should be added to the RHS, but since they are zero there is nothing to add. However, one thing ...


4

You appear (and correct me if I am wrong, i do not wish to presume) to be confused on the notion of what is an order of convergence. The order of convergence is the order at which your approximate solution (in this case, your simulation) approaches the theoretical solution. In your case, the theoretical order of convergence of the centered scheme is 2. ...


4

Your proposed discretisation appears to be consistent, but wouldn't normally be interpreted as a finite volume discretisation. Indeed, it looks a lot more like a finite difference method with a slightly unusual stencil. Further, this looks like a nonlinear hyperbolic, rather than elliptic PDE, which fits with your desire to use FVM. Generically, you can ...


3

Using finite element methods, you also have the choice of discretizing the weak formulation directly, using two types of approaches: using H²-conforming elements like the Argyris element or the Hsieh-Clough-Tocher on triangles or the Bogner-Fox-Schmidt element on quadrilaterals. Using regular H¹-conforming elements in conjunction with the so ...


3

I think that one key point to understand the answers is that, with the parabolic PDE that you wrote, we have some control on the "quantity of $u$", i.e. $\int_{\Omega} u $: If you integrate the PDE, you have $$\int_{\Omega} u_t = \int_{\Omega} \Delta u + \int_{\Omega} f(x,t) = \int_{\Omega} f(x,t)$$ because $\int_{\Omega} \Delta u = \int_{\partial \Omega} \...


3

Yes, with Dirichlet boundary conditions you always have exponential convergence to the staedy state. Any PDE book will have a proof. For a nice explanation from a numerical perspective, see chapter 2 of LeVeque's FDM book.


3

NCLab has three (simple) FEM modules - linear elasticity, electrostatics, and general linear second-order PDE with constant coefficients https://nclab.com/pde. Higher-order elements are enabled but not hp-adaptivity yet. You may still encounter bugs.


3

Since you are creating the mesh specifically for this right-hand side in order to have a correct solution, it would make more sense to have the Dirac located on a node instead: On the one hand, this is likely to be easier to realize in the mesh generator, and on the other hand the assembly is easier: Since you are using continuous piecewise linear elements ...


3

Rewrite your problem as $$\begin{cases} \Delta u(x)=\frac{f(x)}{k_1}& \text{ in } \Omega_1\\ \Delta u(x)=\frac{f(x)}{k_2}& \text{ in } \Omega_2\\ u=0 & \text{ in } \Gamma \end{cases}$$ you have to add a condition on an interface, for example a jump condition, e.g. $$u^1(x)-u^2(x)=\alpha(x)\quad\text{and}\quad u^1_n(x)-u^2_n(y)=\beta(x) ,\...


3

Your statement suggests that you assume that because you have a term $\nabla\cdot k \nabla T$ in the PDE, that the derivatives $\nabla T$ need to be continuous. But that's not true. You take derivatives not of $\nabla T$ but of $k\nabla T$ (i.e., the heat flux), and this quantity is continuous [1]. So it is differentiable almost everywhere, and that's all ...


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