<< Click to Display Table of Contents >> Navigation:  Technical Notes > Avoid Discontinuities

Discontinuities can cause serious numerical difficulty.  This is most glaringly true in time-dependent problems, but can be a factor in steady-state problems as well.

 

Steady-State

 

The finite element model used in FlexPDE assumes that all variables are continuous throughout the problem domain.  This follows from the fact that the mesh nodes that sample the values of the variables are shared between the cells that they adjoin.  Internally, the solution variables are interpolated by low-order polynomials over each cell of the finite element mesh.  A discontinuous change in boundary conditions along the boundary path, particularly between differing VALUE conditions, will require intense mesh refinement to resolve the transition.  

 

Whenever possible, use RAMP, URAMP, SWAGE, part of a sine or supergaussian, or some other smooth function to make a transition in value conditions over a physically meaningful distance.

 

If the quantity you have chosen as a system variable is in fact expected to be discontinuous at an interface, consider choosing a different variable which is continuous, and from which the real variable can be computed.

 

Time-Dependent

 

It is a common tendency in posing problems for numerical solution to specify initial conditions or boundary conditions as discontinuous functions, such as "at time=2 seconds, the boundary temperature is raised instantaneously to 200 degrees."  A little thought will reveal that such statements are totally artificial.  They violate the constraints of physics, and they pose impossible conditions for numerical solution.  Not quite so obvious is the case where a boundary condition is applied at the start of the problem which is inconsistent with the initial values.  This is in fact a statement that "at time=0 the boundary temperature is raised instantaneously to a new value", and so is the same as the statement above.

 

To raise a temperature "instantaneously" requires an infinite heat flux.  To move a material position "instantaneously" requires an infinite force.  In the real world, nothing happens "instantaneously".  Viscosity diffuses velocity gradients, elastic deformation softens displacement velocities, thermal diffusion smoothes temperature changes.  At some scale, all changes in nature are smooth.

 

In the mathematical view, the Fourier transform of a step function is (1/frequency).  This means that a discontinuity excites an infinite spectrum of spatial and temporal frequencies, with weights that diminish quite slowly at higher frequencies.  An "accurate" numerical model of such a system would require an infinite number of nodes and infinitesimal time steps, to satisfy sampling requirements of two samples per cycle.  Any frequency components for which the sampling requirement is not met will be modeled wrong, and will cause oscillations or inaccuracies in the solution.

 

How then have numerical solutions been achieved to these problems over the decades?  The answer is that artificial numerical diffusion processes have secretly filtered the frequency spectrum of the solution to include only low-frequency components.  Or the answers have been wrong.  Right enough to satisfy the user, and wrong enough to satisfy the calculation.  

 

It is useful in this context to note that the effect of a diffusion term D*div(grad(U)) is to apply an attenuation of 1/(1+D*K*K) to the K-th frequency component of U.  Conversely, any side effect of a numerical approximation which damps high frequency components is similar to a diffusion operator in the PDE.

 

We have attempted in FlexPDE to eliminate as many sources of artificial solution behavior as possible.  Automatic timestep control and adaptive gridding are mechanisms which try to follow accurately the solution of the posed PDE.  Discontinuities cannot be accurately modeled, and are therefore, strictly speaking, ill-posed problems.  They cause tiny timesteps and intense mesh refinement in the early phases, causing long running times.

 

What can be done?  

•Start your problem with initial conditions which are self-consistent; this means the values should correspond to a steady state solution with some set of boundary conditions. If you cannot by inspection determine these values, use the INITIAL EQUATIONS facility or a steady-state FlexPDE run with TRANSFER to precompute the initial values.  See initialeq.pde and smoothing_discontinuities.pde in the Samples|Usage|Sequenced_Equations folder.

•Use RAMP, URAMP, SWAGE or other smooth function of time to turn the source value on over a meaningful interval of time.

•Whenever possible, instead of an instantaneously applied value condition, use a flux boundary condition which reflects the maximum physical initial flux that could arise from such a step condition (see the sample problem SAMPLES|Applications|Misc|Diffusion.pde for an example).

•Volume source functions and Natural boundary conditions are not as sensitive as direct conditions on the variables, because they appear in the numerical solution as integrals over some interval, and are thus somewhat "smoothed".

 

It may seem like an imposition that we should require such adulteration of your pure PDE, but the alternative is that we apply these adulterations behind your back, in unknown quantities and with unknown affect on your solution.  At least this way, you're in control.