Tag Info

30

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}\... 11 Curved boundaries are covered in most CFD books, e.g., Chapter 11 of Wesseling or Chapter 8 of Ferziger and Peric. While not a fundamental theoretical problem, the practical complexity of implementing boundary conditions for high-order methods on curved boundaries is a significant reason for interest in more geometrically-flexible methods such as the finite ... 11 In deal.II (http://www.dealii.org -- disclaimer: I'm one of the principal authors of that library), we do eliminate whole rows and columns, and it is not too expensive overall. The trick is to use the fact that the sparsity pattern is typically symmetric, so you know which rows you need to look into when eliminating a whole column. The better approach, in ... 10 The best way to do this is (as you said) to just use the definition of periodic boundary conditions and set up your equations correctly from the start using the fact that u(0)=u(1). In fact, even more strongly, periodic boundary conditions identify x=0 with x=1. For this reason, you should only have one of these points in your solution domain. An open ... 10 The problem you describe, how to prescribe non-reflecting or absorbing boundary conditions when solving partial differential equations (PDE) has been extensively studied. For complex (e.g. nonlinear) systems of PDE and general boundaries, it can be quite challenging and is something of an ongoing research problem. You can find many references on the topic, ... 10 I have found that I must keep the value of dt near dx or the results become unstable This behavior you have noticed is known as the Courant–Friedrichs–Lewy (CFL) condition. Indeed, there are timestep restrictions that are related to spatial refinement and choice of time integration method. Due to the diffusion term in 13a, the timestep must satisfy an ... 9 The solution of the equation you are looking for is in the space H^1 of functions that have one weak derivative, but in 2d and 3d, this does not imply that the solution is in fact continuous. As a consequence, it is not possible to define the value of a solution at individual points, and equally consequentially, it is mathematically not possible to ... 8 Answering your last question first, do people actually use FDM for curved boundary nowadays I'd say the answer is no. In the commercial CFD world, 2nd order accurate finite volume schemes are the de-facto industry standard. One of the advantages of FV (and finite element/discontinuous galerkin approaches Jed mentioned) over FD is the much more natural ... 8 Kirill's answer is a good method but it is hampered that the fact that you need to tune it to a specific range in energy, so it is not appropriate if you have a broad-band wavepacket you need to absorb. An alternative approach is exterior complex scaling, which is reviewed well in Infinite-range exterior complex scaling as a perfect absorber in time-... 8 For the sake of simplicity let me slightly change your notation. Let u, v, w be the components in the x, y, z directions of a true (polar) vector, and y=0 the symmetry plane. For a true vector, symmetry conditions across the xz plane are: \begin{eqnarray} u(x,y,z) &= u(x, -y, z) \\ -v(x,y,z) &= v(x, -y, z) \\ w(x,y,z) &= w(x,-y,z) \end{... 8 It is standard procedure for general case like that. You can enforce condition for electrode by Lagrange multipliers, L(u,\lambda) = \int_\Omega \sigma u_{,i} u_{,i} \, \textrm{d}V - \int_{\partial\Omega_\sigma} u n_i j_i \textrm{d}S + \lambda\left(\int_E \sigma u_{,n}\, \textrm{d}S - I \right) $$that for discretised problem is$$ \left[ \begin{array}{...

7

Typically, you would decompose your solution into $u = u_D + u_0$, where $u_D$ satisfies inhomogeneous Dirichlet conditions. You can then solve for $u_0 \in H^1_0$ subject to $$\Delta u_0 = f-\Delta u_D.$$ You would then recover a variational formulation over $H^1_0$ again. Note that $u_D$ is non-unique. To remedy this in proofs, $u_D$ is often taken to be ...

7

Such problems (sometimes called lateral Cauchy problems) are in general not well-posed (meaning they either lack a solution, or there are infinitely many of them, or the solution is unstable under perturbations of the boundary conditions). For parabolic (or dissipative) equations, it makes sense to study the stationary limit (simply omit the term $u_t$ in ...

7

There is no shortcut. Just like there is no shortcut to becoming an "engineering expert in a short time". The thing is that to be an expert in civil engineering, you need to understand load analysis, use cases of buildings and bridges, materials, designs, and regulatory issues. For computational science, you need to understand the mathematical background, ...

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

I don't have electrostatics training but your second picture looks OK to me. Your conductive inclusion is essentially shorting out the plates u1 and u2. I was taught that there is no electric field in conductive bodies, and therefore electrostatic potential is constant there. If the conductive body touches both plates u1 and u2 which are forced to ...

6

Let me give a part of your answer, I would need some more indications from your side to answer you fully. So please read, and write some comments so that I can complete my answer. About notations in numerical methods There are a few mistakes in the way you write your equations. These are only details, but for someone who is used to it, it can a bit ...

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

There are multiple ways of implementing Dirichlet boundary conditions when using the finite element method. For each unknown $i$ of the system belonging to the Dirichlet boundary, you can zero out row $i$ of the Jacobian and entry $i$ of the right-hand side. At that rate, you can reformulate the problem for just the unknowns in the interior, which reduces ...

6

$K u$ equals the internal forces only in the linear case. The tangent stiffness matrix, $K$, in a nonlinear problem is normally used in a Newton-Raphson algorithm to calculate updates to the displacement vector as follows: $$K \Delta u = f - f_{internal}$$ $$u_{i+1} = u_i + \Delta u$$ The vector of internal forces, $f_{internal}$ must be calculated from ...

6

From Ablowitz and Zeppetella we know that the analytical solution reads: $$u(x,t)=\frac{1}{1+e^{{(-\frac{5}{6} t+ \frac{\sqrt{6}x}{6}})^2}}$$ Usually, analytical solutions require the imposition of boundary conditions. Are there any used in the derivation of this expression? If so, you must use the same boundary conditions ...

6

The problem with ABC in the form of derivative operators at the boundary is that the exact ABCs are non-local in space and often in time, which makes them hard to use in numerical simulations. You can derive approximate ABCS based on series expansions, but this can be tedious and it may be difficult to get stable schemes. This approach was made famous by a ...

6

Typically, you would add "guard cells", that is (for u) u(-1) and u(n+1) with your notation. Before each integration step: u(n+1) = u(0) u(-1) = u(n) and similarly for the other variables. If you use higher order derivatives, you could also define u(-2) and u(n+2), etc.

6

If $\nabla u$ is of more interest than $u$ itself, it is reasonable to set $\mathbb v := -\sigma \, \nabla u$ and convert your Poisson problem $$-\nabla \cdot (\sigma \, \nabla u) = f \quad \text{in }\Omega, \\ u = g_E \quad \text{on }\partial\Omega_E, \\ \sigma \, \nabla u \cdot \hat{\mathbb n} = g_N \quad \text{on }\partial\Omega_N$$ to a mixed one: $$\... 6 Consider the fact that your linear system can be broken up into two parts: the equations that correspond to the known Dirichlet dofs, and the equations that correspond to the unknown dofs. For the sake of convenience, let's say that your unknown vector is broken up so the first k values of T are known, and the remainder are unknown. Let a "u" subscript ... 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 ... 6 In short, no. Neumann boundary conditions should be specified in terms of traction. This is clear when you move from a strong-form statement of the boundary value problem to a weak-form. Strong form: \begin{align} div(\sigma(\mathbf{u})) + \mathbf{b} &= \mathbf{0} &&\forall \mathbf{x} \in \Omega \\ \mathbf{u} &= \hat{\mathbf{u}}(\mathbf{x}) &... 5 It seems that you are referring to the fluxes that occur at the edge of your domain. Since these are typically known quantities enforced by your boundary conditions, there is no need to solve for those values. There are two common ways to deal with this: Include an equation for each (known) boundary quantity in the form f_i = F(x_i) for x_i\in\partial\... 5 The best thing to do here is to actually reduce the size of your problem after you have assembled the full system matrices. Let's assume you've gotten to that point, so you have \mathbf{\tilde{K}d}=\mathbf{\tilde{F}}  where $\mathbf{K}$ is a $N_{dof} \times N_{dof}$ matrix of stiffness coefficients, $\mathbf{d}$ is your full vector of nodal degrees of ...

5

You could augment your system of ODEs to include one more equation. If you let \begin{align} I(x) = \int_{a}^{x}(A(s) + B(s) + C(s))\,\mathrm{d}s, \end{align} then $I(a) = 0$, $I(b) = N$, $\dot{I}(x) = A(x) + B(x) + C(x)$, and you have another boundary value problem that you can solve in MATLAB using bvp4c.

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