Wolfgang Bangerth schrieb am Dienstag, 23. April 2024 um 16:45:36 UTC+2: On 4/23/24 07:59, Andreas Müsing wrote:
> b) How do I implement this in deal.ii? for (cell=...) for (q=...) for (i=...) for (j=...) local_matrix(i,j) += fe_values.shape_grad(i,q) * fe_values.shape_grad(j,q) * fe_values.quadrature_point(q)[0] * // r fe_values.JxW(q); this was the first approach that I tried, but the results were not what I expected. My answer to your reply on c) will clarify it. > c) When I look at the first term in the PDE's strong form, I can carry out the > first term r derivative by applying the derivative product rule, which results > in terms having a first and a second derivative in r, respectively. Therefore, > the solution should be in any case a linear combination of first derivative > and second derivative values. Can you elaborate? It's true that if you multiply everything out, the equation is of the form d/dr ... + d^2/dr^2 ... = ... but I don't understand how you infer something about the solution then? I have a magnetostatic model, where I calculate the magnetic vector potential A(r,z) in cylinder coordinates, assuming a cylinder symmetric current density J(r,z). For a delta peak excitation in J(r,z), the magnetic field B = curl A is known analytically on the symmetry axis (B(r=0, z)). This delta excitation (in FEM parlance: the rhs vector has a single non-zero component) is a current circle filament in 3D space. For the same current circle filament, I created a separate verification solution via the Biot-Savart law over the whole domain, to get an alternative solution B(r,z). Your implementation suggestion (which I tried before), introduces a factor r in the matrix. This by itself did not lead to the expected spatial solution of A(r,z) and B(r,z), when I compared with B_z(r=0) or the Biot-Savart solution. Then my reasoning was, that the involved integrals contain the expression ∂/ ∂r (r ∂A / ∂r) = ∂A/∂r + r ∂^2/∂r^2 A , and I assumed that the Laplace equation in cylinder coordinates requires an additional first derivative term. > I used tutorial step-63 with success as a starting point, where the > advection-diffusion equation has both a second derivative (diffusion) term, > and a first derivative (=advection) term. Further, I used an analytic known > solution of my problem to compare the results, and I get reasonable solutions > for e = 1 and b = 2. But why beta = 2 and not beta = 1? What is the form of the advection equation you are trying to solve, in cylindrical coordinates? There was a misunderstanding: I am still solving a Laplace problem, but the transition from Cartesian coordinates to cylinder coordinates added a term to the equation which can be *interpreted* as an advection term. The structure of the PDE in step-63 <https://www.dealii.org/current/doxygen/deal.II/step_63.html> is therefore exactly what I need. I played with the parameters epsilon and beta, and for the mentioned epsilon=1.0 and beta = 2.0, the solution matched exactly my analytic results. I am happy with my implementation except of the important question that bothers me: why beta = 2? Maybe I did somewhere else a silly mistake, but most probably I just don't have yet enough insight into FEM and its math details. Thanks for your help! -- The deal.II project is located at http://www.dealii.org/ For mailing list/forum options, see https://groups.google.com/d/forum/dealii?hl=en --- You received this message because you are subscribed to the Google Groups "deal.II User Group" group. To unsubscribe from this group and stop receiving emails from it, send an email to dealii+unsubscr...@googlegroups.com. To view this discussion on the web visit https://groups.google.com/d/msgid/dealii/32682253-16de-4d2d-862a-6883f3813a0fn%40googlegroups.com.
