# Tag Info

15

I think that one of your problems is that (as you observed in your comments) Neumann conditions are not the conditions you are looking for, in the sense that they do not imply the conservation of your quantity. To find the correct condition, rewrite your PDE as $$\frac{\partial \phi}{\partial t} = \frac{\partial}{\partial x}\left( D\frac{\partial \phi}{\... 15 This is a well-framed question and a very useful thing to understand. Korrok is correct to refer you to von Neumann analysis and LeVeque's book. I can add a bit more to that. I'd like to write a detailed answer, but at the moment I only have time for a short one: With \alpha=\beta=1/2, you get a method that is absolutely stable for arbitrarily large ... 10 Generally speaking, you'll want to use an implicit method for parabolic equations (the diffusion part) -- explicit schemes for parabolic PDE need to have a very short timestep to be stable. Conversely, for the hyperbolic part (advection) you'll want an explicit method as it's cheaper and doesn't disrupt the symmetry of the linear system you have to solve by ... 8 Solving partial differential equations with explicit timestepping methods relies on meeting a certain CFL condition for stability. Since you are using a forward Euler timestepping scheme, you must ensure that your timestep is small enough compared to the grid spacing. In the case of the heat equation (Fick's 2nd law is a heat equation), your CFL number C ... 8 It's easy to derive that equation from Fick's law. You have this diffusion equation as:$$\frac{\partial C}{\partial t} = D \nabla^{2} C$$The mean square displacement weighted by the concentration profile is defined as:$$\langle r^{2}(t) \rangle = \int_{\Omega} |\vec{r}|^{2} C(\vec{r},t) d^{3} \vec{r}$$The time-evolution of this mean square ... 5 I think it might be ill-posed, since the time-dependent parts are linearly dependent. If you add your two time-dependent equations together, you get a time-independent equation: (\alpha(x)u_x)_x + (\beta(x) v_x)_x = 0. Do your initial conditions satisfy this equation? 5 This is possible (see [1]) but uncommon, as it requires Monte Carlo moves that alter the current conformations by a very small perturbation. In that setting of "small" Metropolis MC moves, it is usually easier (both in theory and in practice) to just use Molecular Dynamics instead. [1] Kikuchi, K., Yoshida, M., Maekawa, T., & Watanabe, H. (1991). ... 5 It is in some cases possible to map the dynamics obtained in MC simulations to other (more realistic) dynamics, especially for the case of dense colloidal suspensions. The following two papers talk about the problems and caveats of performing such a mapping: http://aip.scitation.org/doi/10.1063/1.3414827 (spherical particles) http://aip.scitation.org/doi/10.... 5 First off, the PDE can be rewritten instead as$$\frac{\partial C}{\partial t} = \frac{\partial}{\partial x}C\frac{\partial C}{\partial x}$$or, by applying the product rule in reverse again, as$$\frac{\partial C}{\partial t} = \frac{1}{2}\frac{\partial^2}{\partial x^2}C^2.$$This equation is often referred to as the porous medium equation or the slow ... 5 Yes, the problem with mixed boundary conditions is well posed. What's not clear to me is this: Why do you approximate the derivative via the two-sides approximation? Shouldn't it be enough to just the following?$$ \frac{f_6 - f_{6-\Delta x}}{\Delta x} = 0. $$In your animations, what is the size of \Delta x? Your curve looks very smooth, which can be ... 4 That's not a simple question in my opinion. It is like if you were doing an experiment in a box and you just want to simulate the center of it. But obviously, the solution in the center of the box depends on what happens in the rest of the box. That is probably why it is a difficult problem to state. I have no definite answer to this question, but I still ... 4 The answer is yes. Drift diffusion = advection-diffusion. You will find this rather frequently: different disciplines, trying to describe the phenomena in their field, assign names to the equations they derive that have meaning for this discipline. Other disciplines come up with the same equations, in different contexts, and name them differently. It's still ... 4 That's a modelling problem, not a CS one. This will depend on what is diffusing and on the reasons why 2D is relevant. You may want to give physical details and ask this on physics.SE. Alternatively if you just want to get your hand in with solving diffusion, then the numbers don't have to reflect a physical fact (especially if you have diffusion only, the ... 4 More generally, the equation you are solving is $$\frac{\partial u}{\partial t} = \nabla \cdot \left( D(x) \nabla u \right) \quad ,$$ where the first nabla operator represents the divergence operator, and the second represents the gradient operator. Your equation is the special case for a one-dimensional cartesian coordinate ... 3 The Crank-Nicolson discretization of this equation will read$$ \frac{T^n-T^{n-1}}{\Delta t} = \frac 12 \left[ \partial_x \left((T^n)^{5/2} \partial_x T^n\right) + \partial_x \left((T^{n-1})^{5/2} \partial_x T^{n-1}\right) \right] $$which is a nonlinear, time-independent, elliptic partial differential equation in T^n. The way to solve such ... 3 You should rearrange the terms so that all of the n+1 terms are together on one side of the equals sign and all of the n terms are on them other. Then you will have a system of non-linear equations like:$$ A(u^{n+1}_i) u^{n+1}_i=b $$Where A is a non-linear matrix. Then you can use something like Newton-Raphson to linearize the system and your ... 3 In statistics the widely used test for checking if the distribution is gaussian is the Jarque-Bera test. Koizumi [1] presents an equivalent test for the multivariate case. I don't know if there something ready to use in C++, Fortran, or Python. However, at first sight, the test looks easy to implement. [1] - Koizumi, Kazuyuki, Naoya Okamoto, and Takashi Seo.... 3 The best method I know of for this problem is by Griffiths. 3 You just add the diffusion along the other dimensions. This superposition from orthogonal directions makes some sense, as they are independent. So, going by wikipedia for Fick's second law of diffusion in 1D:$$ \frac{\partial \psi}{\partial t} = D \frac{\partial^2 \psi}{\partial x^2} $$We extend it to 2d as:$$ \frac{\partial \psi}{\partial t} = D \...

3

There is no reason to believe that two random fields with the same arithmetic mean would yield solutions that have anything to do with each other. In fact, for the case you consider, one might imagine that maybe the harmonic mean is actually a better indicator, but even that is unclear -- it could also be the geometric mean. Apart from this, you have to ...

2

You should never use explicit method for the diffusion equation. Implicit is unconditionally stable and just as easy to implement. Also if you use an implicit method (like backward Euler or Crank-Nicolson) it will not matter how small d is. In fact you could use a dirac delta function if you wanted. As far as speed goes, doing a implicit method with a ...

2

For the heat equation, the wave speed is infinite and this is reflected by the fact that you are solving a globally coupled differential equation if you are using an implicit time stepping scheme.

2

You need to unfold the positions after the simulation is done. I use the following subroutine: subroutine unfold_positions(L, X, Xu) ! Unwinds periodic positions. It is using the positions of particles from the ! first time step as a reference and then tracks them as they evolve possibly ! outside of this box. The X and Xu arrays are of the type X(:, j, i), ...

2

Your first test should be to compute the mean value and covariance matrix of your point sample. If these converge to the correct values as you increase the number of samples, you are, from a practical perspective on the safe side that your points indeed come from the correct distribution. Of course, in practice there are many distributions that have the ...

2

I would do it exactly how you describe, with the two following remarks. I would solve the two evolutions simultaneously. In the way you describe, you will need to store and retrieve all the time steps for the solution $u_A$. Instead, you could solve sequentially for $u_A$ and $u_B$ at each time step, since $u_B$ just requires two solutions for $u_A$. If you ...

2

At first glance, it looks like you are using the method of lines with forward Euler time steps. What WolfgangBangerth is getting at with his example is that even for a simple heat equation, the stability limit of forward Euler (namely, that $|\lambda\Delta{t}| < 1$) combined with the eigenvalues induced by a finite difference approximation of the ...

2

Interesting question. I too would like to know if this is mentioned anywhere explicitly. My guess is maybe this is a little too hard and not too useful to be used/taught. Perhaps someone can give a proper reference. My calculations (see code below) suggest, if I am not wrong, that for the equation  u_t = \mathcal{L}u, \qquad \mathcal{L}u = \nabla\cdot(D\...

2

One of the best ways to test a PDE solver is to use the method of manufactured solutions. Essentially, you modify the PDE (and discretization) by adding a source term that yields an exact solution known in advance. You can then compare your numerical solution against the exact solution for debugging purposes. You should test your numerical solution against a ...

2

You are applying Richardson extrapolation to increasingly-accurate independent sample paths of an OU process. How would this work even for independent identically normally distributed random variables $y_1,y_2,y_3,y_4\sim\mathrm{N}(\mu,\sigma^2)$? These have already "converged" to the true distribution $\mathrm{N}(\mu,\sigma^2)$, and you might expect the ...

2

If you want to obtain the diffusion coefficient for 2D from 3D you have to take a look how the 2D diffusion equation is derived from 3D equation. It will be integrated with respect to one of the spatial dimensions. Diffusion coefficient is a $3 \times 3$ tensor in 3D, in 2D it will be $2 \times 2$ tensor. So in my opinion you cannot simply take the value ...

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