Пресс-релиз популярных книг

Авторы: 111 А Б В Г Д Е Ж З И Й К Л М Н О П Р С Т У Ф Х Ц Ч Ш Щ Э Ю Я

Книги: 164 А Б В Г Д Е Ж З И Й К Л М Н О П Р С Т У Ф Х Ц Ч Ш Щ Э Ю Я

На сайте 111 авторов, 92 книг, 72 статей, 5913 глав.

Parabolic Problems

Back

___ Introduction

The _nite element method may be used to solve time_dependent problems as well as

steady ones_ This e_ort involves both parabolic and hyperbolic partial di_erential sys_

tems_ Problems of parabolic type involve di_usion and dissipation while hyperbolic

problems are characterized by conservation of energy and wave propagation_ Simple

one_dimensional heat conduction and wave propagation equations will serve as model

problems of each type_

Example ______ The one_dimensional heat conduction equation

ut _ puxx_ a _ x _ b_ t _ __ ___ _ a

where p is a positive constant called the di_usivity_ is of parabolic type_ Initial_boundary

value problems consist of determining u_x_ t satisfying ___ _ a given the initial data

u_x_ _ _ u__x_ a _ x _ b_ ___ _ b

and appropriate boundary data_ e_g__

pux_a_ t _ __u_a_ t _ ___t_ pux_b_ t _ __u_b_ t _ ___t_ ___ _ c

As with elliptic problems_ boundary conditions without the pux term are called Dirichlet

conditions those with _i _ __ i _ __ _ are Neumann conditions and those with both

terms present are called Robin conditions_ The problem domain is open in the time

direction t thus_ unlike elliptic systems_ this problem is evolutionary and computation

continues in t for as long as there is interest in the solution_

Example ______ The one_dimensional wave equation

utt _ c_uxx_ a _ x _ b_ t _ __ ___ __a

_ Parabolic Problems

where c is a constant called the wave speed_ is a hyperbolic partial di_erential equation_

Initial_boundary value problems consist of determining u_x_ t satisfying ___ __a given

the initial data

u_x_ _ _ u__x_ ut_x_ _ _ u_ __x_ a _ x _ b_ ___ __b

and boundary data of the form ___ _ c_ Small transverse vibrations of a taut string

satisfy the wave equation_ In this case_ u_x_ t is the transverse displacement of the

string and c_ _ T___ T being the applied tension and _ being the density of the string_

We_ll study parabolic problems in this chapter and hyperbolic problems in the next_

We shall see that there are two basic _nite element approaches to solving time_dependent

problems_ The _rst_ called the method of lines_ uses _nite elements in space and ordinary

di_erential equations software in time_ The second uses _nite element methods in both

space and time_ We_ll examine the method of lines approach _rst and then tackle space_

time _nite element methods_

___ Semi_Discrete Galerkin Problems_ The Method

of Lines

Let us consider a parabolic problem of the form

ut _ L_u_ _ f_x_ y_ _x_ y _ __ t _ __ _____ a

where L is a second_order elliptic operator_ In two dimensions_ u would be a function of

x_ y_ and t and L_u_ could be the very familiar

L_u_ _ __puxx _ _puyy _ qu_ _____ b

Appropriate initial and boundary conditions would also be needed_ e_g__

u_x_ y_ _ _ u__x_ y_ _x_ y _ _ _ __ _____ c

u_x_ y_ t _ _x_ y_ t_ _x_ y _ _E_ _____ d

pun _ _u _ __ _x_ y _ _N_ _____ e

We construct a Galerkin formulation of _____ in space in the usual manner thus_ we

multiply _____ a by a suitable test function v and integrate the result over _ to obtain

_v_ ut _ _v_L_u_ _ _v_ f_

____ Semi_Discrete Galerkin Problems _

As usual_ we apply the divergence theorem to the second_derivative terms in L to reduce

the continuity requirements on u_ When L has the form of _____ b_ the Galerkin problem

consists of determining u _ H_

E _ _t _ _ such that

_v_ ut _ A_v_ u _ _v_ f_ _ v__ _ _u __ _v _ H_

_ _ t _ __ ______a

The L_ inner product_ strain energy_ and boundary inner product are_ as with elliptic

problems_

_v_ f _ ZZ _

vfdxdy_ ______b

A_v_ u _ ZZ _

_p_vxux _ vyuy _ vqu_dxdy_ ______c

and

_ v_pun __ Z__N

vpunds_ ______d

The natural boundary condition _____ e has been used to replace pun in the boundary

inner product_ Except for the presence of the _v_ ut term_ the formulation appears to

the same as for an elliptic problem_

Initial conditions for ______a are usually determined by projection of the initial data

_____ c either in L_

_v_ u _ _v_ u__ _v _ H_

_ _ t _ __ ______a

or in strain energy

A_v_ u _ A_v_ u__ _v _ H_

_ _ t _ __ ______b

Example ______ We analyze the one_dimensional heat conduction problem

ut _ _puxx _ f_x_ t_ _ _ x _ _ t _ __

u_x_ _ _ u__x_ _ _ x _ _

u___ t _ u_ _ t _ __ t _ __

thoroughly in the spirit that we did in Chapter for a two_point boundary value problem_

A Galerkin form of this heat_conduction problem consists of determining u _ H_

satisfying

_v_ ut _ A_v_ u _ _v_ f_ _v _ H_

_ _ t _ __

_ Parabolic Problems

0 = x x x x

U(x,t)

0 1 j N-1 N x = 1

N-1

Figure ____ _ Mesh for the _nite element solution of Example ____ _

_v_ u _ _v_ u__ _v _ H_

_ _ t _ __

where

A_v_ u _ Z _

vxpuxdx_

Boundary terms of the form ______d disappear because v _ _ at x _ __ with Dirichlet

data_

We introduce a mesh on _ _ x _ as shown in Figure ____ and choose an approxi_

mation U of u in a _nite_dimensional subspace SN

_ of H_

_ having the form

U_x_ t _

N__

Xj__

cj_t_j_x_

Unlike steady problems_ the coe_cients cj_ j _ _ __ _ _ _ _N_ _ depend on t_ The Galerkin

_nite element problem is to determine U _ SN

_ such that

__j_ Ut _ A__j_ U _ __j_ f_ t _ __

__j_ U _ __j_ u__ t _ __ j _ _ __ _ _ _ _N _ _

Let us chose a piecewise_linear polynomial basis

_k_x _ __

x_xk__

xk_xk__

_ if xk__ _ x _ xk

xk___x

xk___xk

_ if xk _ x _ xk__

__ otherwise

This problem is very similar to the one_dimensional elliptic problem considered in Section

___ so we_ll skip several steps and also construct the discrete equations by vertices rather

than by elements_

____ Semi_Discrete Galerkin Problems _

Since _j has support on the two elements containing node j we have

A__j_ U _ Z xj

xj__

jpUxdx _ Z xj__

jpUxdx

where _ _ _ d_ _dx_ Substituting for _j and Ux

A__j_ U _ Z xj

xj__

p_x_

cj _ cj__

dx _ Z xj__

xj _

hj__

p_x_

cj__ _ cj

hj__

where

hj _ xj _ xj___

Using the midpoint rule to evaluate the integrals_ we have

A__j_ U _

pj____

_cj _ cj__ _

pj____

hj__

_cj__ _ cj

where pj____ _ p_xj_____

Similarly_

__j_ Ut _ Z xj

xj__

_jUtdx _ Z xj__

_jUtdx

__j_ Ut _ Z xj

xj__

_j_ _ cj___j__ _ _ cj_jdx _ Z xj__

_j_ _ cj_j _ _ cj___j__dx

where __ _ d_ _dt_ Since the integrands are quadratic functions of x they may be

integrated exactly using Simpson_s rule to yield

__j_ Ut _

_ _ cj__ _ __ cj _

hj__

__ _ cj _ _ cj___

Finally_

__j_ f _ Z xj

xj__

_jf_xdx _ Z xj__

_jf_xdx_

Although integration of order one would do_ we_ll_ once again_ use Simpson_s rule to

obtain

__j_ f _

__fj____ _ fj _

hj__

_fj _ _fj_____

We could replace fj____ by the average of fj__ and fj to obtain a similar formula to the

one obtained for __j_ Ut thus_

__j_ f _

_fj__ _ _fj _

hj__

__fj _ fj___

Combining these results yields the discrete _nite element system

_ _ cj__ _ __ cj _

hj__

__ _ cj _ _ cj__ _

pj____

_cj _ cj__ _

pj____

hj _ __

_cj__ _ cj

_ Parabolic Problems

_fj__ _ _fj _

hj__

__fj _ fj___ j _ _ __ _ _ _ _N _ _

_We have dropped the _ and written the equation as an equality_

If p is constant and the mesh spacing h is uniform_ we obtain

_ _ cj__ _ __ cj _ _ cj__ _

_cj__ _ _cj _ cj__ _

_fj__ _ _fj _ fj___

j _ _ __ _ _ _ _N _ _

The discrete systems may be written in matrix form and_ for simplicity_ we_ll do so for

the constant coe_cient_ uniform mesh case to obtain

M_c _ Kc _ l ______a

where

M _

_______

_ _ _

_ ______b

K _

_______

_ _

_ _ _

_ _

_ ______c

l _

_____

f_ _ _f_ _ f_

___

fN__ _ _fN__ _ fN

_ ______d

c _ _c__ c__ _ _ _ _ cN___T _ ______e

The matricesM_ K_ and l are the global mass matrix_ the global sti_ness matrix_ and

the global load vector_ Actually_ M has little to do with mass and should more correctly

be called a global dissipation matrix however_ we_ll stay with our prior terminology_

In practical problems_ element_by_element assembly should be used to construct global

matrices and vectors and not the nodal approach used here_

The discrete _nite element system ______ is an implicit system of ordinary di_erential

equations for _c_ The mass matrix M can be _lumped_ by a variety of tricks to yield an

____ Semi_Discrete Galerkin Problems _

explicit ordinary di_erential system_ One such trick is to approximate __j_ Ut by using

the right_rectangular rule on each element to obtain

__j_ Ut _ Z xj

xj__

_j_ _ cj___j__ _ _ cj_jdx _ Z xj__

_j_ _ cj_j _ _ cj___j__dx _ hcj_

The resulting _nite element system would be

hI _c _ Kc _ l_

Recall _cf_ Section ____ that a one_point quadrature rule is satisfactory for the conver_

gence of a piecewise_linear polynomial _nite element solution_

With the initial data determined by L_ projection onto SN

E _ we have

__j_ U___ _ _ __j_ u__ j _ _ __ _ _ _ _N _ _

Numerical integration will typically be needed to evaluate __j_ u_ and we_ll approximate

it in the manner used for the loading term __j_ f_ Thus_ with uniform spacing_ we have

Mc__ _ u_ _

_____

_ _u_

_ u_

_ _u_

_ u_

___

N__ _ _u_

N__ _ u_

_ ______f

If the initial data is consistent with the trivial Dirichlet boundary data_ i_e__ if u_ _ H_

then the above system reduces to

cj__ _ u__xj_ j _ _ __ __ _ _ _ _N _ _

Had we solved the wave equation ___ __ instead of the heat equation ___ _ using a

piecewise_linear _nite element basis_ we would have found the discrete system

M_c _ Kc _ _ ______

with p in ______c replaced by c__

The resulting initial value problems _IVPs for the ordinary di_erential equations

_ODEs ______a or ______ typically have to be integrated numerically_ There are several

excellent software packages for solving IVPs for ODEs_ When such ODE software is used

with a _nite element or _nite di_erence spatial discretization_ the resulting procedure is

called the method of lines_ Thus_ the nodes of the _nite elements appear to be _lines_

in the time direction and_ as shown in Figure _____ for a one_dimensional problem_ the

temporal integration proceeds along these lines_

_ Parabolic Problems

0 = x0 x1 xj xN-1 N x = 1

Figure ______ _Lines_ for a method of lines integration of a one_dimensional problem_

Using the ODE software_ solutions are calculated in a series of time steps ___ t___

_t__ t___ _ _ _ _ Methods fall into two types_ Those that only require knowledge of the so_

lution at time tn in order to obtain a solution at time tn__ are called one_step methods_

Correspondingly_ methods that require information about the solution at tn and several

times prior to tn are called multistep methods_ Excellent texts on the subject are available

___ __ __ ___ One_step methods are Runge_Kutta methods while the common multistep

methods are Adams or backward di_erence methods_ Software based on these methods

automatically adjusts the time steps and may also automatically vary the order of accu_

racy of a class of methods in order to satisfy a prescribed local error tolerance_ minimize

computational cost_ and maintain numerical e_ciency_

The choice of a one_step or multistep method will depend on several factors_ Gener_

ally_ Runge_Kutta methods are preferred when time integration is simple relative to the

spatial solution_ Multistep methods become more e_cient for complex nonlinear prob_

lems_ Implicit Runge_Kutta methods may be e_cient for problems with high_frequency

oscillations_ The ODEs that arise from the _nite element discretization of parabolic

problems are _sti__ ___ __ so backward di_erence methods are the preferred multistep

methods_

Most ODE software ___ __ __ addresses _rst_order IVPs of the explicit form

y_ _t _ f _t_ y_t_ y__ _ y__ ______

Second_order systems such as ______ would have to be written as a _rst_order system by_

e_g__ letting

d _ _c

____ Semi_Discrete Galerkin Problems _

and_ hence_ obtaining

_ _c

M_d _ _ _ d

_Kc __

Unfortunately_ systems having the form of ______a or the one above are implicit and

would require inverting or lumping M in order to put them into the standard explicit

form _______ Inverting M is not terribly di_cult when M is constant or independent

of t however_ it would be ine_cient for nonlinear problems and impossible when M is

singular_ The latter case can occur when_ e_g__ a heat conduction and a potential problem

are solved simultaneously_

Codes for di_erential_algebraic equations _DAEs directly address the solution of im_

plicit systems of the form

f _t_ y_t_ y_ _t _ __ y__ _ y__ ______

One of the best of these is the code DASSL written by Petzold ____ DASSL uses variable_

step_ variable_order backward di_erence methods to solve problems without needingM__

to exist_

Let us illustrate these concepts by applying some simple one_step schemes to problems

having the forms _____ or _______ However_ implementation of these simple methods

is only justi_ed in certain special circumstances_ In most cases_ it is far better to use

existing ODE software in a method of lines framework_

For simplicity_ we_ll assume that all boundary data is homogeneous so that the bound_

ary inner product in ______a vanishes_ Selecting a _nite_dimensional space SN

_ H_

__ we

then determine U as the solution of

_V_Ut _ A_V_U _ _V_ f_ _v _ SN

_ _ ______

Evaluation leads to ODEs having the form of ______a regardless of whether or not the

system is one_dimensional or the coe_cients are constant_ The actual matrices M and K

and load vector l will_ of course_ di_er from those of Example ____ in these cases_ The

systems ______a or ______ are called semi_discrete Galerkin equations because time has

not yet been discretized_

We discretize time into a sequence of time slices _tn_ tn___ of duration _t with tn _

n_t_ n _ __ _ _ _ _ _ For this discussion_ no generality is lost by considering uniform time

steps_ Let_

u_x_ tn be the exact solution of the Galerkin problem ______a at t _ tn_

U_x_ tn be the exact solution of the semi_discrete Galerkin problem ______ at t _ tn_

Un_x be the approximation of U_x_ tn obtained by ODE software_

_ Parabolic Problems

cj_tn be the Galerkin coe_cient at t _ tn thus_ for a one_dimensional problem

U_x_ tn _

N__

Xj__

cj_tn_j_x_

For a Lagrangian basis_ cj_tn _ U_xj_ tn_

cnj

be the approximation of cj_tn obtained by ODE software_ For a one_dimensional

problem

Un_x _

N__

Xj__

cnj

_j_x_

We suppose that all solutions are known at time tn and that we seek to determine

them at time tn___ The simplest numerical scheme for doing this is the forward Euler

method where ______ is evaluated at time tn and

Ut_x_ tn _

Un___x _ Un_x

_ ______

A simple Taylor_s series argument reveals that the local discretization error of such an

approximation is O__t_ Substituting ______ into ______ yields

_V_

Un__ _ Un

_ A_V_Un _ _V_ fn_ _v _ SN

_ _ _____ _a

Evaluation of the inner products leads to

cn__ _ cn

_ Kncn _ ln_ _____ _b

We have allowed the sti_ness matrix and load vector to be functions of time_ The mass

matrix would always be independent of time for di_erential equations having the explicit

form of _____ a as long as the spatial _nite element mesh does not vary with time_

The ODEs _____ _a_b are implicit unless M is lumped_ If lumping were used and_ e_g__

M _ hI then cn__ would be determined as

cn__ _ cn _

_ln _ Kncn__

Assuming that cn is known_ we can determine cn__ by inverting M_

Using the backward Euler method_ we evaluate ______ at tn__ and use the approxi_

mation ______ to obtain

_V_

Un__ _ Un

_ A_V_Un__ _ _V_ fn___ _v _ SN

_ _ _____ a

____ Semi_Discrete Galerkin Problems

and

cn__ _ cn

_ Kn__cn__ _ ln___ _____ b

The backward Euler method is implicit regardless of whether or not lumping is used_

Computation of cn__ requires inversion of

M_ Kn___

The most popular of these simple schemes uses a weighted average of the forward and

backward Euler methods with weights of _ _ and __ respectively_ Thus_

_V_

Un__ _ Un

_ _ _ _A_V_Un _ _A_V_Un__ _ _ _ __V_ fn _ __V_ fn___

_V _ SN

_ _ _____ _a

and

cn__ _ cn

_ _ _ _Kncn _ _Kn__cn__ _ _ _ _ln _ _ln___ _____ _b

The forward and backward Euler methods are recovered by setting _ _ _ and _ respec_

tively_

Let us regroup terms involving cn and cn__ in _____ _b to obtain

_M_ __tKn___cn__ _ _M_ _ _ __tKn_cn __t__ _ _ln _ _ln____ _____ _c

Thus_ determination of cn__ requires inversion of

M_ __tKn___

In one dimension_ this system would typically be tridiagonal as with Example ____ _ In

higher dimensions it would be sparse_ Thus_ explicit inversion would never be performed_

We would just solve the sparse system _____ _c for cn___

Taylor_s series calculations reveal that the global discretization error is

kc_tn _ cnk _ O__t

for almost all choices of _ _ ___ _ ____ The special choice _ _ __ yields the Crank_Nicolson

method which has a discretization error

kc_tn _ cnk _ O__t__

The foregoing discussion involved one_step methods_ Multistep methods are also used

to solve time_dependent _nite element problems and we_ll describe them for an ODE in

_ Parabolic Problems

the implicit form _______ The popular backward di_erence formulas _BDFs approximate

y_t in ______ by a k th degree polynomial Y_t that interpolates y at the k _ times

tn___i_ i _ __ _ _ _ _ _ k_ The derivative _ y is approximated by _Y

_ The Newton backward

di_erence form of the interpolating is most frequently used to represent Y ___ ___ but

since we_re more familiar with Lagrangian interpolation we_ll write

y_t _ Y_t _

Xi__

yn___iNi_t_ t _ _tn___k_ tn____ _____ _a

where

Ni_t _

Y j___j__i

t _ tn___j

tn___i _ tn___j

_ _____ _b

The basis _____ _b is represented by the usual Lagrangian shape functions _cf_ Section

____ so Ni_tn___j _ ij _

Assuming yn___i_ i _ _ __ _ _ _ _ k_ to be known_ the unknown yn__ is determined by

collocation at tn___ Thus_ using ______

f _tn___Y_tn___ _Y

_tn__ _ __ _____ _

Example ______ The simplest BDF formula is obtained by setting k _ in _____ _ to

obtain

Y_t _ yn__N__t _ ynN__t_

N__t _

t _ tn

tn__ _ tn

_ N__t _

t _ tn__

tn _ tn__

Di_erentiating Y_t

_t _

yn__ _ yn

tn__ _ tn

thus_ the numerical method _____ _ is the backward Euler method

f _tn___ yn___

yn__ _ yn

tn__ _ tn

_ __

Example ______ The second_order BDF follows by setting k _ _ in _____ _ to get

Y_t _ yn__N__t _ ynN__t _ yn__N__t

N__t _

_t _ tn_t _ tn__

__t_ _ N__t _

_t _ tn___t _ tn__

__t_ _

N__t _

_t _ tn___t _ tn

__t_ _

where time steps are of duration _t_

____ Finite Element Methods in Time _

Di_erentiating and setting t _ tn__

__tn__ _

__t

_ _N

__tn__ _ _

_ _N

__tn__ _

__t

Thus_

_tn__ _

_yn__ _ _yn _ yn__

__t

and the second_order BDF is

f _tn___ yn___

_yn__ _ _yn _ yn__

__t

_ __

Applying this method to ______a yields

_cn__ _ _cn _ cn__

__t

_ Kn__cn__ _ ln___

Thus_ computation of cn__ requires inversion of

__t

_ K_

Backward di_erence formulas through order six are available ___ __ __ __ ___

___ Finite Element Methods in Time

It is_ of course_ possible to use the _nite element method in time_ This can be done

on space_time triangular or quadrilateral elements for problems in one space dimension

on hexahedra_ tetrahedra_ and prisms in two space dimensions and on four_dimensional

parallelepipeds and prisms in three space dimensions_ However_ for simplicity_ we_ll focus

on the time aspects of the space_time _nite element method by assuming that the spatial

discretization has already been performed_ Thus_ we_ll consider an ODE system in the

form ______a and construct a Galerkin problem in time by multiplying it by a test

function w _ L_ and integrating on _tn_ tn___ to obtain

_w_M_cn _ _w_Kcn _ _w_ ln_ _w _ L__tn_ tn____ _____ a

where the L_ inner product in time is

_w_ cn _ Z tn__

wT cdt_ _____ b

Only _rst derivatives are involved in ______a thus_ neither the trial space for c nor the

test space for w have to be continuous_ For our initial method_ let us assume that c_t

is continuous at tn_ By assumption_ c_tn is known in this case and_ hence_ w_tn _ __

_ Parabolic Problems

Example ______ Let us examine the method that results when c_t and w_t are linear

on _tn_ tn____ We represent c_t in the manner used for a spatial basis as

c__ _ cnNn__ _ cn__Nn____ ______a

where

Nn__ _

_ _

_ Nn____ _

_ _

______b

are hat functions in time and

_ _

_t _ tn _ tn__

______c

de_nes the canonical element in time_ The test function

w _ Nn____ _ _ _ _ _ _ _ _T ______d

vanishes at tn __ _ _ and is linear on _tn_ tn___

Transforming the integrals in _____ a to __ _ using ______c and using ______a_b_d

yields

_w_M_cn _

_ Z _

_ _

cn__ _ cn

d__

_w_Kcn _

_ Z _

_ _

K_cn _ _

_ cn__ _ _

_d__

_Again_ we have written equality instead of _ for simplicity_ Assuming that M and K

are independent of time_ we have

_w_M_cn _ M

cn__ _ cn

_w_Kcn _

K_cn _ _cn___

Substituting these into _____ a

cn__ _ cn

K_cn _ _cn__ _

_ Z _

_ _

l__ d_ ______a

or_ if l is approximated like c_

cn__ _ cn

K_cn _ _cn__ _

_ln _ _ln___ ______b

Regrouping terms

_M_

_tK_cn__ _ _M_

_tK_cn _

_t_ln _ _ln____ ______c

____ Finite Element Methods in Time _

we see that the piecewise_linear Galerkin method in time is a weighted average scheme

_____ _c with _ _ ____ Thus_ at least to this low order_ there is not much di_erence be_

tween _nite di_erence and _nite element methods_ Other similarities appear in Problem

at the end of this section_

Low_order schemes such as _____ _ are popular in _nite element packages_ Our pref_

erence is for BDF or implicit Runge_Kutta software that control accuracy through au_

tomatic time step and order variation_ Implicit Runge_Kutta methods may be derived

as _nite element methods by using the Galerkin method _____ with higher_order trial

and test functions_ Of the many possibilities_ we_ll examine a class of methods where the

trial function c_t is discontinuous_

Example ______ Suppose that c_t is a polynomial on _tn_ tn___ with jump disconti_

nuities at tn_ n _ __ When we need to distinguish left and right limits_ we_ll use the

notation

cn_ _ lim

___

c_tn _ __ cn_ _ lim

___

c_tn _ __ ______a

With jumps at tn_ we_ll have to be more precise about the temporal inner product _____ b

and we_ll de_ne

_u_ vn_ _ lim

___ Z tn____

tn__

uvdt_ _u_ vn_ _ lim

___ Z tn____

tn__

uvdt_ ______b

The inner product _u_ vn_ may be a_ected by discontinuities in functions at tn_ but

_u_ vn_ only involves integrals of smooth functions_ In particular_

_u_ vn_ _ _u_ vn_ when u_t and v_t are either continuous or have jump discon_

tinuities at tn

_u_ vn_ exists and _u_ vn_ _ _ when either u or v are proportional to the delta

function _t _ tn and

_u_ vn_ doesn_t exist while _v_ un_ _ _ when both u and v are proportional to

_t _ tn_

Suppose_ for example_ that v_t is continuous at tn and u_t _ _t _ tn_ Then

_u_ vn_ _ lim

___ Z tn____

tn__

_t _ tnv_tdt _ v_tn_

The delta function can be approximated by a smooth function that depends on _ as was

done in Section ___ to help explain this result_

Let us assume that w_t is continuous and write c_t in the form

c_t _ cn_ _ _ c_t _ cn__H_t _ tn ______a

_ Parabolic Problems

where

H_t _ _ if t _ _

__ otherwise

______b

is the Heaviside function and c is a polynomial in t_

Di_erentiating

_c_t _ _ c_t _ cn___t _ tn _ _

c_tH_t _ tn_ ______c

With the interpretation that inner products in _____ are of type _______ assume that

w_t is continuous and use ______ in _____ a to obtain

wT _tnM_tn_cn_ _ cn_ _ _w_M_

cn_ _ _w_K cn_ _ _w_ ln__ _w _ H__ ______

The simplest discontinuous Galerkin method uses a piecewise constant _p _ _ basis

in time_ Such approximations are obtained from ______a by selecting

c_t _ cn_ _ c_n__ __

Testing against the constant function

w_t _ _ _ _ _ _ _ _ _T

and assuming that M and K are independent of t_ ______ becomes

M_c_n__ _ _ cn_ _ Kc_n__ __t _ Z tn__

l_tdt_

The result is almost the same as the backward Euler formula _____ b except that the

load vector l is averaged over the time step instead of being evaluated at tn___

With a linear _p _ approximation for c_t_ we have

c_t _ cn_Nn_t _ c_n__ _Nn___t

where Nn_i_ i _ __ _ are given by ______b_ Selecting the basis for the test space as

wi_t _ Nn_i_t_ _ _ _ _ _ _ _T _ i _ __ _

assuming thatMand K are independent of t_ and substituting the above approximations

into _______ we obtain

M_cn_ _ cn_ _

M_c_n__ _ _ cn_ _

K__cn_ _ c_n__ _ _ Z tn__

Nnl_tdt

____ Finite Element Methods in Time _

and

M_c_n__ _ _ cn_ _

K_cn_ _ _c_n__ _ _ Z tn__

Nn__l_tdt_

Simplifying the expressions and assuming that l_t can be approximated by a linear

function on _tn_ tn__ yields the system

cn_ _ c_n__ _

_ _ cn_ _

K__cn_ _ c_n__ _ _

__ln _ l_n__ __

c_n__ _ _ cn_

K_cn_ _ _c_n__ _ _

_ln _ _l_n__ __

This pair of equations must be solved simultaneously for the two unknown solution vectors

cn_ and c_n__ __ This is an implicit Runge_Kutta method_

Problems

_ Consider the Galerkin method in time with a continuous basis as represented by

_____ _ Assume that the solution c_t is approximated by the linear function

______a_c on _tn_ tn__ as in Example ____ _ but do not assume that the test space

w_t is linear in time_

_ _ Specifying

w__ _ ___ _ _ _ _ _ _ _ _T

and assuming that M and K are independent ot t_ show that _____ a is the

weighted average scheme

_M_ __tK_cn__ _ _M_ _ _ __tK_cn __t__ _ _ln _ _ln___

with

_ _ R_

__ ___ N_

n____ d_

__ ___ d_

When di_erent trial and test spaces are used_ the Galerkin method is called a

Petrov_Galerkin method_

___ The entire e_ect of the test function __t is isolated in the weighting factor __

Furthermore_ no integration by parts was performed_ so that __t need not be

continuous_ Show that the choices of __t listed in Table ____ correspond to

the cited methods_

__ The discontinuous Galerkin method may be derived by simultaneously discretizing

a partial di_erential system in space and time on _ _ _t _ n__ t_n__ __ This form

may have advantages when solving problems with rapid dynamics since the mesh

may be either moved or regenerated without concern for maintaining continuity

_ Parabolic Problems

Scheme _ _

Forward Euler _____ _b _ _ _ _

Crank_Nicolson _____ _b __ !_

Crank_Nicolson _____ _b !_

Backward Euler _____ b _ _ _

Galerkin ______ N_

n____ _!_

Table ____ _ Test functions _ and corresponding methods for the _nite element solution

of ______a with a linear trial function_

between time steps_ Using ______a as a model spatial _nite element formulation_

assume that test functions v_x_ y_ t are continuous but that trial functions u_x_ y_ t

have jump discontinuities at tn_ Assume Dirichlet boundary data and show that

the space_time discontinuous Galerkin form of the problem is

_v_ utST _ _v___ tn_ u___ tn_ _ u___ tn_ _ AST _v_ u _ _v_ fST _

_v _ H_

___ _ _tn__ t_n__ __

where

_v_ uST _ Z t_n____

tn_ ZZ _

vudxdydt

and

AST _v_ u _ _vx_ puxST _ _vy_ puyST _ _v_ quST _

In this form_ the _nite element problem is solved on the three_dimensional strips

_ _ _tn__ t_n__ __ n _ __ _ _ _ _ _

___ Convergence and Stability

In this section_ we will study some theoretical properties of the discrete methods that

were introduced in Sections ___ and ____ Every _nite di_erence or _nite element scheme

for time integration should have three properties_

_ Consistency_ the discrete system should be a good approximation of the di_erential

equation_

__ Convergence_ the solution of the discrete system should be a good approximation

of the solution of the di_erential equation_

__ Stability_ the solution of the discrete system should not be sensitive to small per_

turbations in the data_

____ Convergence and Stability _

Somewhat because they are open ended_ _nite di_erence or _nite element approxi_

mations in time can be sensitive to small errors_ e_g__ introduced by round o__ Let us

illustrate the phenomena for the weighted average scheme _____ _c

_M_ __tK_cn__ _ _M_ _ _ __tK_cn __t__ _ _ln _ _ln____ _____

We have assumed_ for simplicity_ that K and M are independent of time_

Sensitivity to small perturbations implies a lack of stability as expressed by the fol_

lowing de_nition_

De_nition ______ A _nite di_erence scheme is stable if a perturbation of size kk in_

troduced at time tn remains bounded for subsequent times t _ T and all time steps

_t _ _t__

We may assume_ without loss of generality_ that the perturbation is introduced at

time t _ __ Indeed_ it is common to neglect perturbations in the coe_cients and con_ne

the analysis to perturbations in the initial data_ Thus_ in using De_nition ____ _ we

consider the solution of the related problem

_M_ __tK_"cn__ _ _M_ _ _ __tK_"cn __t__ _ _ln _ _ln____

"c_ _ c_ _ __

Subtracting _____ from the perturbed system

_M_ __tK__n__ _ _M_ _ _ __tK__n_ __ _ __ ______a

where

_n _ "cn _ cn_ ______b

Thus_ for linear problems_ it su_ces to apply De_nition ____ to a homogeneous version

of the di_erence scheme having the perturbation as its initial condition_ With these

restrictions_ we may de_ne stability in a more explicit form_

De_nition ______ A linear di_erence scheme is stable if there exists a constant C _ _

which is independent of _t and such that

k_nk _ Ck__k ______

as n__ _t _ __ t _ T_

__ Parabolic Problems

Both De_nitions ____ and _____ permit the initial perturbation to grow_ but only

by a bounded amount_ Restricting the growth to _nite times t _ T ensures that the

de_nitions apply when the solution of the di_erence scheme cn _ as n__ When

applying De_nition ______ we may visualize a series of computations performed to time

T with an increasing number of time steps M of shorter_and_shorter duration _t such

that T _ M_t_ As _t is decreased_ the perturbations n_ n _ _ __ _ _ _ _M_ should settle

down and eventually not grow to more than C times the initial perturbation_

Solutions of continuous systems are often stable in the sense that c_t is bounded for

all t _ __ In this case_ we need a stronger de_nition of stability for the discrete system_

De_nition ____ _ The linear di_erence scheme _____ is absolutely stable if

k_nk _ k__k_ ______

Thus_ perturbations are not permitted to grow at all_

Stability analyses of linear constant coe_cient di_erence equations such as ______

involve assuming a perturbation of the form

_n _ __nr_ ______

Substituting into ______a yields

_M_ __tK___n__r _ _M_ _ _ __tK___nr_

Assuming that _ __ _ and M_ __tK is not singular_ we see that _ is an eigenvalue and

r is an eigenvector of

_M_ __tK____M_ _ _ __tK_rk _ _krk_ k _ _ __ _ _ _ _N_ ______

Thus_ _n will have the form ______ with _ _ _k and r _ rk when the initial perturbation

__ _ rk_ More generally_ the solution of ______a is the linear combination

_n _

Xk__

k__knrk ______a

when the initial perturbation has the form

__ _

Xk__

krk_ ______b

Using ______a_ we see that ______ will be absolutely stable when

j_kj _ _ k _ _ __ _ _ _ _N_ ______

____ Convergence and Stability _

The eigenvalues and eigenvectors of many tridiagonal matrices are known_ Thus_ the

analysis is often straight forward for one_dimensional problems_ Analyses of two_ and

three_dimensional problems are more di_cult however_ eigenvalue_eigenvector pairs are

known for simple problems on simple regions_

Example ______ Consider the eigenvalue problem ______ and rearrange terms to get

_M_ __tK__krk _ _M_ _ _ __tK_rk

__k _ Mrk _ ___k_ _ _ _ ___tKrk

_Krk _ _kMrk

where

_k _

__k_ _ _ _ ___t

Thus_ _k is an eigenvalue and rk is an eigenvector of _M__K_

Let us suppose that M and K correspond to the mass and sti_ness matrices of the

one_dimensional heat conduction problem of Example ____ _ Then_ using ______b_c_ we

have

_____

_ _

_ _ _

_ _

_____

rk_

___

rk_N__

_kh

_____

_ _ _

_____

rk_

___

rk_N__

The di_usivity p and mesh spacing h have been assumed constant_ Also_ with Dirichlet

boundary conditions_ the dimension of this system is N _ rather than N_

It is di_cult to see in the above form_ but writing this eigenvalue_eigenvector problem

in component form

_rj__ _ _rj _ rj__ _

_kh

_rj__ _ _rj _ rj___ j _ _ __ _ _ _ _N _ _

we may infer that the components of the eigenvector are

rkj _ sin

k_j

_ j _ _ __ _ _ _ _N _ _

This guess of rk may be justi_ed by the similarity of the discrete eigenvalue problem to

a continuous one however_ we will not attempt to do this_ Assuming it to be correct_ we

substitute rkj into the eigenvalue problem to _nd

_sin

k__j _

N _ _ sin

k_j

_ sin

k__j _

__ Parabolic Problems

_kh

_sin

k__j _

_ _sin

k_j

_ sin

k__j _

_ j _ _ __ _ _ _ _N _ _

But

sin

k__j _

_ sin

k__j _

_ _ sin

k_j

cos

and

_cos

N _ sin

k_j

_kh

_cos

_ _ sin

k_j

Hence_

_k _ __p

h___cos k__N _

cos k__N _ ___

With cos k__N ranging on __ _ __ we see that _ _p_h_ _ _k _ __ Determining _k in

terms of _k

_k _

_ _k_ _ __t

_ _k__t

_ _

_k_t

_ _k__t

We would like j_kj _ for absolute stability_ With _k _ __ we see that the requirement

that _k _ is automatically satis_ed_ Demanding the _k _ _ yields

j_kj_t_ _ __ _ __

If _ _ __ then _ __ _ _ and the above inequality is satis_ed for all choices of _k and

_t_ Methods of this class are unconditionally absolutely stable_ When _ _ ___ we have

to satisfy the condition

p_t

h_ _

__ _ __

If we view this last relation as a restriction of the time step _t_ we see that the forward

Euler method __ _ _ has the smallest time step_ Since all other methods listed in Table

____ are unconditionally stable_ there would be little value in using the forward Euler

method without lumping the mass matrix_ With lumping_ the stability restriction of the

forward Euler method actually improves slightly to p_t_h_ _ ___

Let us now turn to a more general examination of stability and convergence_ Let_s

again focus on our model problem_ determine u _ H_

_ satisfying

_v_ ut _ A_v_ u _ _v_ f_ _v _ H_

_ _ t _ __ ______a

_v_ u _ _v_ u__ _v _ H_

_ _ t _ __ ______b

The semi_discrete approximation consists of determining U _ SN

_ H_

_ such that

_V_Ut _ A_V_U _ _V_ f_ _V _ SN

_ _ t _ __ _____ _a

____ Convergence and Stability __

_V_U _ _V_ u__ _V _ SN

_ _ t _ __ _____ _b

Trivial Dirichlet boundary data_ again_ simpli_es the analysis_

Our _rst result establishes the absolute stability of the _nite element solution of the

semi_discrete problem _____ _ in the L_ norm_

Theorem ______ Let _ SN

_ satisfy

_V_ t _ A_V_ _ __ _V _ SN

_ _ t _ __ _____ a

_V_ _ _V_ __ _V _ SN

_ _ t _ __ _____ b

Then

k___ __ tk_ _ k_k__ t _ __ _____ c

Remark __ With _x_ t being the di_erence between two solutions of _____ _a satis_

fying initial conditions that di_er by __x_ the loading _V_ f vanishes upon subtraction

_as with _______

Proof_ Replace V in _____ a by to obtain

__ t _ A__ _ __

dtkk_

_ A__ _ __

Integrating

k___ __ tk_

_ k___ __ _k_

_ _ _ Z t

A__ d__

The result _____ c follows by using the initial data _____ b and the non_negativity of

A__ _

We_ve discussed stability at some length_ so now let us turn to the concept of conver_

gence_ Convergence analyses for semi_discrete Galerkin approximations parallels the lines

of those for elliptic systems_ Let us_ as an example_ establish convergence for piecewise_

linear solutions of _____ _ to solutions of _______

Theorem ______ Let SN

_ consist of continuous piecewise_linear polynomials on a family

of uniform meshes _h characterized by their maximum element size h_ Then there exists

a constant C _ _ such that

max

t____T ku _ Uk_ _ C_ _ j log

h_ jh_ max

t____T kuk__ _____ _

__ Parabolic Problems

Proof_ Create the auxiliary problem_ determine W _ SN

_ such that

__V_W_ ___ __ _ _ A_V_W___ __ _ _ __ _V _ SN

_ _ _ _ ___ t_ _____ _a

W_x_ y_ t _ E_x_ y_ t _ U_x_ y_ t _ " U_x_ y_ t_ _____ _b

where "U _ SN

_ satis_es

A_V_ u___ __ _ _ " U___ __ _ _ __ _V _ SN

_ _ _ _ ___ T__ _____ _c

We see that W satis_es a terminal value problem on _ _ _ _ t ant that " U satis_es an

elliptic problem with _ as a parameter_

Consider the identity

_W_ E _ _W__ E _ _W_ E_ _

Integrate and use _____ _b

kE___ __ tk_

_ _W_ E___ __ _ _ Z t

__W__ E _ _W_ E_ _d__

Use _____ _a with V replaced by E

kE___ __ tk_

_ _W_ E___ __ _ _ Z t

_A_W_ E _ _W_ E_ _d__ _____ _

Setting v in ______ and V in _____ _ to W and subtracting yields

_W_ u_ _ U_ _ A_W_ u _ U _ __ _ _ __

_W_ u _ U__ _ __ _ _ __

Add these results to _____ _ and use _____ _b to obtain

kE___ __ tk_

_ _W_ ____ __ _ _ Z t

_A_W_ _ _ _W_ __ _d__

where

_ _ u _ " U_

The _rst term in the integrand vanishes by virtue of _____ _c_ The second term is

integrated by parts to obtain

kE___ __ tk_

_ _W_ ____ __ t _ Z t

_W__ _d__ _____ _a

____ Convergence and Stability __

This result can be simpli_ed slightly by use of Cauchy_s inequality _j_W_ V j _ kWk_kV k_

to obtain

kE___ __ tk_

_ _ kW___ __ tk_k____ __ tk_ _ Z t

_ kW_k_k_k_d__ _____ _b

Introduce a basis on SN

_ and write W in the standard form

W_x_ y_ _ _

Xj__

cj__ _j_x_ y_ _____ _

Substituting _____ _ into _____ _a and following the steps introduced in Section ____ we

are led to

_M_c _ Kc _ __ _____ _a

where

Mij _ __i_ _j_ _____ _b

Kij _ A__i_ _j_ i_ j _ _ __ _ _ _ _N_ _____ _c

Assuming that the sti_ness matrix K is independent of _ _ _____ _a may be solved exactly

to show that _cf_ Lemmas ____ and _____ which follow

kW___ __ _k_ _ kE___ __ tk__ _ _ _ _ t_ _____ _a

Z t

_ kW_k_d_ _ C_ _ j log

h_ jkE___ __ tk__ _____ _b

Equation _____ _a is used in conjunction with _____ _b to obtain

kE___ __ tk_

_ _ _kE___ __ tk_ _ Z t

_ kW_k_d_ max

_____t k____ __ _k__

Now_ using _____ _b

kE___ __ tk_ _ C_ _ j log

h_ j max

_____t k____ __ _k__ _____ _

Writing

u _ U _ u _ "U

_ " U _ U _ _ _ E

and taking an L_ norm

ku _ Uk_ _ k_k_ _ kEk__

__ Parabolic Problems

Using _____ _

ku _ Uk_ _ C_ _ j log

h_ j max

_____t k____ __ _k__ _______a

Finally_ since _ satis_es the elliptic problem _____ _c_ we can use Theorem _____ to

write

k____ __ _k_ _ Ch_ku___ __ _k__ _______b

Combining _______a and _______b yields the desired result _____ __

The two results that were used without proof within Theorem _____ are stated as

Lemmas_

Lemma ______ Under the conditions of Theorem ______ there exists a constant C _ _

such that

A_V_ V _

h_ kV k_

_ _V _ SN

_ _ ______

Proof_ The result can be inferred from Example ____ however_ a more formal proof is

given by Johnson ____ Chapter __

Instead of establishing _____ _b_ we_ll examine a slightly more general situation_ Let

c be the solution of

M_c _ Kc _ __ t _ __ c__ _ c__ _______

The mass and sti_ness matrices M and K are positive de_nite_ so we can diagonalize

________ In particular_ let be a diagonal matrix containing the eigenvalues of M__K

and R be a matrix whose columns are the eigenvectors of the same matrix_ i_e__

M__KR _ R_ _______a

Further let

d_t _ R__c_t_ _______b

Then _______ can be written in the diagonal form

_ d _ _ _______a

by multiplying it by _MR__ and using _______a_b_ The initial conditions generally

remain coupled through _______a_b_ i_e__

d__ _ d_ _ R__c__ _______b

With these preliminaries_ we state the desired result_

____ Convection_Di_usion Systems __

Lemma ______ If d_t is the solution of ______ then

j_d j _ jdj _

Cjd_j

_ t _ __ _______a

where jdj _ pdTd_ If_ in addition_

max

____

j_j

j_j _

h_ _______b

then

Z T

_j_d j _ jdjdt _ C_ _ j log

h_ jjd_j_ _______c

Proof_ cf_ Problem _

Problems

_ Prove Lemma ______

__ Convection_Diusion Systems

Problems involving convection and di_usion arise in #uid #ow and heat transfer_ Let us

consider the model problem

ut _ _ _ ru _ r _ _pru _____ a

where _ _ ____ ___T is a velocity vector_ Written is scalar form_ _____ a is

ut _ __ux _ __uy _ _puxx _ _puyy_ _____ b

The vorticity transport equation of #uid mechanics has the form of _____ _ In this case_

u would represent the vorticity of a two_dimensional #ow_

If the magnitude of _ is small relative to the magnitude of the di_usivity p_x_ y_

then the standard methods that we have been studying work _ne_ This_ however_ is not

the case in many applications and_ as indicated by the following example_ standard _nite

element methods can produce spurious results_

Example _____ _ __ Consider the steady_ one_dimensional_ convection_di_usion equa_

tion

__u__ _ u_ _ __ _ _ x _ _ ______a

with Dirichlet boundary conditions

u__ _ _ u_ _ __ ______b

__ Parabolic Problems

The exact solution of this problem is

u_x _ _

e____x __ _ e____

_ e____ _ ______c

If _ _ _ _ then_ as shown by the solid line in Figure ____ _ the solution features

a boundary layer near x _ _ At points removed from an O__ neighborhood of x _ _

the solution is smooth with u _ _ Within the boundary layer_ the solution rises sharply

from its unit value to u _ _ at x _ _

0 0.2 0.4 0.6 0.8 1

−4

−3

−2

−1

2 N odd

N even

Figure ____ _ Solutions of ______ with _ _ ____ The exact solution is shown as a solid

line_ Piecewise_linear Galerkin solutions with __ and _element meshes are shown as

dashed and dashed_dotted lines_ respectively _ __

The term _u__ is di_usive while the term u_ is convective_ With a small di_usivity

__ convection dominates di_usion outside of the narrow O__ boundary layer_ Within

this layer_ di_usion cannot be neglected and is on an equal footing with convection_

This simple problem will illustrate many of the di_culties that arise when _nite element

methods are applied to convection_di_usion problems while avoiding the algebraic and

geometric complexities of more realistic problems_

Let us divide ___ _ into N elements of width h _ _N_ Since the solution is slowly

varying over most of the domain_ we would like to choose h to be signi_cantly larger than

____ Convection_Di_usion Systems __

the boundary layer thickness_ This could introduce large errors within the boundary layer

which we assume can be reduced by local mesh re_nement_ This strategy is preferable to

the alternative of using a _ne mesh everywhere when the solution is only varying rapidly

within the boundary layer_

Using a piecewise_linear basis_ we write the _nite element solution as

U_x _

Xj__

cj_j_x_ c_ _ _ cN _ __ ______a

where

_k_x ___

x_xk__

xk_xk__

_ if xk__ _ x _ xk

xk___x

xk___xk

_ if xk _ x _ xk__

__ otherwise

_ ______b

The coe_cients c_ and cN are constrained so that U_x satis_es the essential boundary

conditions ______b_

The Galerkin problem for ______ consists of determining U_x _ SN

_ such that

____

i_ U_ _ __i_ U_ _ __ i _ _ __ _ _ _ _N _ _ ______a

Since this problem is similar to Example ____ _ we_ll omit the development and just write

the inner products

___

i_ U_ _

_ci__ _ _ci _ ci___ ______b

__i_ U_ _

ci__ _ ci__

_ ______c

Thus_ the discrete _nite element system is

_ _

ci__ _ _ci _ _ _

ci__ _ __ i _ _ __ _ _ _ _N _ _ ______d

The solution of this second_order_ constant_coe_cient di_erence equation is

ci _ _

_ _i

_ _N _ i _ __ _ _ _ _ _N_ ______e

_ _

_ h___

_ ______f

The quantity h___ is called the cell Peclet or cell Reynolds number_ If h___ _ _ then

_ _ _

_ O__

_ _ eh__ _ O__

__ Parabolic Problems

which is the correct solution_ However_ if h___ _ _ then _ _ _ and

ci _ _ if i is even

__ if i is odd

when N is odd and

ci _ _N _ i_N_ if i is even

O_ ___ if i is odd

when N is even_ These two obviously incorrect solutions are shown with the correct

results in Figure ____ _

Let us try to remedy the situation_ For simplicity_ we_ll stick with an ordinary di_er_

ential equation and consider a two_point boundary value problem of the form

L_u_ _ __u__ _ _u_ _ qu _ f_ _ _ x _ _ ______a

u__ _ u_ _ __ ______b

Let us assume that u_ v _ H_

_ with u_ and v_ being continuous except_ possibly_ at

x _ _ _ ___ _ Multiplying ______a by v and integrating the second derivative terms by

parts yields

_v_L_u_ _ A_v_ u _ __u_v_x__ ______a

where

A_v_ u _ __v__ u_ _ _v_ _u_ _ _v_ qu_ ______b

_Q_x__ _ lim

___

_Q__ _ _ Q__ _ __ ______c

We must be careful because the _strain energy_ A_v_ u is not an inner product since

A_u_ u need not be positive de_nite_ We_ll use the inner product notation here for

convenience_

Integrating the _rst two terms of ______b by parts

_v_L_u_ _ _L__v__ u _ ___v_u _ u_v _ _vu_x__

or_ since u and v are continuous

_v_L_u_ _ _L__v__ u _ ___v_u _ u_v_x__ ______a

The di_erential equation

L__v_ _ __v__ _ __v_ _ qv_ ______b

with the boundary conditions v__ _ v_ _ _ is called the adjoint problem and the

operator L__ _ is called the adjoint operator_

____ Convection_Di_usion Systems _

De_nition ______ A Green_s function G___ x for the operator L_ _ is the continuous

function that satis_es

L__G___ x_ _ __Gxx _ __Gx _ qG _ __ x _ ___ _ _ ___ _ ______a

G___ _ _ G___ _ _ ______b

_Gx___ x_x__ _ _

_ ______c

Evaluating ______a with v_x _ G___ x while using ______a_ _____ and assuming that

u__x _ H____ gives the familiar relationship

u__ _ _L_u__G___ _ _ Z _

G___ xf_xdx_ ______a

A more useful expression for our present purposes is obtained by combining ______a and

______a with v_x _ G___ x to obtain

u__ _ A_u_G___ __ ______b

As usual_ Galerkin and _nite element Galerkin problems for ______a would consist of

determining u _ H_

_ or U _ SN

_ H_

_ such that

A_v_ u _ _v_ f_ _v _ H_

_ _ _____ _a

and

A_V_U _ _V_ f_ _v _ SN

_ _ _____ _b

Selecting v _ V in _____ _a and subtracting _____ _b yields

A_V_ e _ __ _v _ SN

_ _ _____ _c

where

e_x _ u_x _ U_x_ _____ _d

Equation ______b did not rely on the continuity of u__x hence_ it also holds when u

is replaced by either U or e_ Replacing u by e in ______b yields

e__ _ A_e_G___ __ _____ a

__ Parabolic Problems

Subtacting _____ _c

e__ _ A_e_G___ _ _ V _ _____ b

Assuming that A_v_ u is continuous in H__ we have

je__j _ Ckek_kG___ _ _ V k__ _____ c

Expressions _____ b_c relate the local error at a point _ to the global error_ Equation

_____ c also explains superconvergence_ From Theorem _____ we know that kek_ _

O_hp when SN consists of piecewise polynomials of degree p and u _ Hp___ The test

function V is also an element of SN however_ G___ x cannot be approximated to the

same precision as u because it may be less smooth_ To elaborate further_ consider

kG___ _ _ V k_

Xj__

kG___ _ _ V k_

where

kuk_

_j _ Z xj

xj__

__u__ _ u__dx_

If _ _ _xk___ xk_ k _ _ __ _ _ _ _N_ then the discontinuity in Gx___ x occurs on some

interval and G___ x cannot be approximated to high order by V _ If_ on the other hand_

_ _ xk_ k _ __ _ _ _ _ _N_ then the discontinuity in Gx___ x is con_ned to the mesh and

G___ x is smooth on every subinterval_ Thus_ in this case_ the Green_s function can be

approximated to O_hp by the test function V and_ using _____ c_ we have

u_xk _ Ch_p_ k _ __ _ _ _ _ _N_ _____ _

The solution at the vertices is converging to a much higher order than it is globally_

Equation _____ c suggests that there are two ways of minimizing the pointwise error_

The _rst is to have U be a good approximation of u and the second is to have V be a

good approximation of G___ x_ If the problem is not singularly perturbed_ then the two

conditions are the same_ However_ when _ _ _ the behavior of the Green_s function is

hardly polynomial_ Let us consider two simple examples_

Example _____ ____ Consider ______ in the case when __x _ __ x _ ___ __ Balancing

the _rst two terms in ______a implies that there is a boundary layer near x _ thus_

at points other than the right endpoint_ the small second derivative terms in ______ may

be neglected and the solution is approximately

_u_

R _ quR _ f_ _ _ x _ _ uR__ _ __

____ Convection_Di_usion Systems __

where uR is called the reduced solution_ Near x _ the reduced solution must be

corrected by a boundary layer that brings it from its limiting value of uR_ to zero_

Thus_ for _ _ _ _ _ the solution of ______ is approximately

u_x _ uR_x _ uR_ e____x ___ ___

Similarly_ the Green_s function ______ has boundary layers at x _ _ and x _ ___ At

points other than these_ the second derivative terms in ______a may be neglected and

the Green_s function satis_es the reduced problem

___GR_ _ qGR _ __ x _ ___ _ _ ___ _ GR___ x _ C___ _ GR___ _ __

Boundary layer jumps correct the reduced solution at x _ _ and x _ _ and determine an

asymptotic approximation of G___ x as

G___ x _ c__ GR___ x _ GR___ _e____ x___ if x _ _

e__x__ ___ ___ if x _ _

The function c__ is given in Flaherty and Mathon ____

Knowing the Green_s function_ we can construct test functions that approximate it

accurately_ To be speci_c_ let us write it as

G___ x _

Xj__

G___ xj_j_x _____ _

where _j_x_ j _ __ _ _ _ _ _N_ is a basis_ Let us consider ______ and ______ with _ _ __

x _ ___ __ Approximating the Green_s function for arbitrary _ is di_cult_ so we_ll restrict

_ to xk_ k _ __ _ _ _ _ _N_ and establish the goal of minimizing the pointwise error of

the solution_ Mapping each subinterval to a canonical element_ the basis _j_x_ x _

_xj___ xj__ is

_j_x _ $ __

x _ xj

_____ _a

where

$ __s ___

__e_____s_

__e__ _ if _ _ s _ _

e__s_e__

__e__ _ if _ _ s _

__ otherwise

_____ _b

where

_ _

h _

_____ _c

__ Parabolic Problems

−1 −0.8 −0.6 −0.4 −0.2 0 0.2 0.4 0.6 0.8 1

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

Figure ______ Canonical basis element $ __s for _ _ __ __ and __ _increasing steepness_

is the cell Peclet number_ The value of _ will remain unde_ned for the moment_ The

canonical basis element $ __s is illustrated in Figure ______ As _ _ _ the basis _____ _b

becomes the usual piecewise_linear hat function

$ __s _

___

_ s_ if _ _ s _ _

_ s_ if _ _ s _

__ otherwise

As ___ _____ _b becomes the piecewise_constant function

$ __s _ _ if _ _ s _ _

__ otherwise

The limits of this function are nonuniform at s _ _ _ __

We_re now in a position to apply the Petrov_Galerkin method with U _ SN

_ and

V _ $ SN

_ to _______ The trial space SN will consist of piecewise linear functions and_ for

the moment_ the test space will remain arbitrary except for the assumptions

_j_x _ H____ __ _j_xk _ jk_ Z _

$ __sds _ _ j_ k _ _ __ _ _ _ _N _ _

_____ _

____ Convection_Di_usion Systems __

The Petrov_Galerkin system for ______ is

____

i_ U_ _ __i_ _U_ _ __i_ qU _ __i_ f_ i _ _ __ _ _ _ _N _ _ _____ _

Let us use node_by_node evaluation of the inner products in _____ __ For simplicity_ we_ll

assume that the mesh is uniform with spacing h and that _ and q are constant_ Then

____

i_ U_ _

h Z _

$ ___s $ U__sds

where $ U_s is the mapping of U_x onto the canonical element _ _ s _ _ With a

piecewise linear basis for $U

and the properties noted in _____ _ for _j _ we _nd

____

i_ U_ _ _

_ci_ _____ _a

We_ve introduced the central di_erence operator

ci _ ci____ _ ci____ _____ _b

for convenience_ Thus_

_ci _ _ci _ ci__ _ _ci _ ci___ _____ _c

Considering the convective term_

___i_ U_ _ _ Z _

$ __s$U

__sds _ ___ _ ____ci _____ _a

where _ is the averaging operator

_ci _ _ci____ _ ci_______ _____ _b

Thus_

_ci _ __ci _ _ci__ _ ci_____ _____ _c

Additionally_

_ _ _Z _

_ $ __s _ $ ___s_ds _____ _d

Similarly

q__i_ U _ qh Z _

$ __s $ U_sds _ qh_ _ __ _ ___ci _____ _a

__ Parabolic Problems

where

_ Z _

__ jsj $ __sds_ _____ _b

_ _ _Z _

s $ __sds_ _____ _c

Finally_ if f_x is approximated by a piecewise_linear polynomial_ we have

__i_ f _ h_ _ __ _ ___fi _______

where fi _ f_xi_

Substituting _____ _a_ _____ _a_ _____ _a_ and _______ into _____ _ gives a di_erence

equation for ck_ k _ _ __ _ _ _ _N _ _ Rather than facing the algebraic complexity_ let us

continue with the simpler problem of Example ____ _

Example ______ Consider the boundary value problem _______ Thus_ q _ f_x _ _ in

_____ ________ and we have

____

i_ U_ _ ___i_ U_ _ _

_ci _ ___ _ ____ci_ i _ _ __ _ _ _ _N _ _ ______ a

or_ using _____ _c_ _____ _c_ and _____ _c

__ _

_ci__ _ _ci _ ci__ _

ci__ _ ci__

_ __ i _ _ __ _ _ _ _N _ _ ______ b

This is to be solved with the boundary conditions

c_ _ _ cN _ __ ______ c

The exact solution of this second_order constant_coe_cient di_erence equation is

ci _ _

_ _i

_ _N _ i _ __ _ _ _ _ _N_ _______a

where

_ _

_ _ ___ _

_ _______b

In order to avoid the spurious oscillations found in Example ____ _ we_ll insist that

_ _ __ Using _______b_ we see that this requires

_ _ sgn_ _

_ _______c

Some speci_c choices of _ follow_

____ Convection_Di_usion Systems __

_ Galerkin_s method_ _ _ __ In this case_

$ __s _ $__s _

_ jsj

Using ________ this method is oscillation free when

j_j

_ _

From _____ _c_ this requires h _ _j___j_ For small values of j___j_ this would be

too restrictive_

__ Il_in_s scheme_ In this case_ $ __s is given by _____ _b and

_ _ coth

_ _

This scheme gives the exact solution at element vertices for all values of __ Either

this result or the use of _______c indicates that the solution will be oscillation free

for all values of __ This choice of _ is shown with the function _ ___ in Figure

______

__ Upwind dierencing_ _ _ sgn__ When _ _ __ the shape function $ __s is the

piecewise constant function

$ __s _ _ if _ _ s _ _

__ otherwise

This function is discontinuous however_ _nite element solutions still converge_

With _ _ _ _______b becomes

_ _

__ _ __

___

In the limit as ___ we have _ _ _ thus_ using _______a

ci _ _ ___N_i _ i _ __ _ _ _ _ _N_ _ _ _

This result is a good asymptotic approximation of the true solution_

Examining ______ as a _nite di_erence equation_ we see that positive values of _ can

be regarded as adding dissipation to the system_

This approach can also be used for variable_coe_cient problems and for nonuniform

mesh spacing_ The cell Peclet number would depend on the local value of _ and the

mesh spacing in this case and could be selected as

_j _

hj _j

_______

__ Parabolic Problems

0 1 2 3 4 5 6 7 8 9 10

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

Figure ______ The upwinding parameter _ _ coth ___ _ ___ for Il_in_s scheme _upper

curve and the function _ ___ _lower curve vs_ __

where hj _ xj _ xj__ and _j is a characteristic value of __x when x _ _xj___ xj_ e_g__

_j _ __j_____ Upwind di_erencing is too di_usive for many applications_ Il_in_s scheme

o_ers advantages_ but it is di_cult to extend to problems other than _______

The Petrov_Galerkin technique has also been applied to transient problems of ther

form _____ however_ the results of applying Il_in_s scheme to transient problems have

more di_usion than when it is applied to steady problems_

Example _____ ____ Consider Burgers_s equation

_uxx _ uux _ __ _ _ x _ _

with the Dirichlet boundary conditions selected so that the exact solution is

u_x _ tanh

_ x

Burgers_s equation is often used as a test problem because it is a nonlinear problem with

a known exact solution that has a behavior found in more complex problems_ Flaherty

___ solved problems with h__ _ __ ___ and N _ __ using upwind di_erencing and Il_in_s

scheme _the Petrov_Galerkin method with the exponential weighting given by _____ _b_

____ Convection_Di_usion Systems __

h__ Maximum Error

Upwind Exponential

_ __ __ ______

___ _______ ____ __

Table ____ _ Maximum pointwise errors for the solution of Example _____ using upwind

di_erencing __ _ sgn_ and exponential weighting __ _ coth ___ _ ___ ____

The cell Peclet number _______ used

_j ___

U_xj_ if _Uj____ _ _

_U_xj_____ if _Uj____ _ _

U_xj _ _ if _Uj____ _ _

The nonlinear solution is obtained by iteration with the values of U_x evaluated at the

beginning of an iterative step_

The results for the pointwise error

jej_ _ max

__j_N ju_xj _ U_xjj

are shown in Table ____ _ The value of h__ _ _ is approximately where the great_

est di_erence between upwind di_erencing __ _ sgn_ and exponential weighting __ _

coth ___ _ ___ exists_ Di_erences between the two methods decrease for larger and

smaller values of h___

The solution of convection_di_usion problems is still an active research area and much

more work is needed_ This is especially the case in two and three dimensions_ Those

interested in additional material may consult Roos et al_ _ ___

Problems

_ Consider ______ when __x __ q_x _ __ x _ ___ _ ____

_ _ Show that the solution of ______ is asymptotically given by

u_x _

f_x

q_x _ uR__e_xpq__ __ _ uR_ e____x pq__ ___

Thus_ the solution has O_p_ boundary layers at both x _ _ and x _ _

___ In a similar manner_ show that the Green_s function is asymptotically given

G___ x _

____q_xq______ _ e____x pq__ ___ if x _ _

e__x__ pq__ ___ if x _ _

The Green_s function is exponentially small away from x _ __ where it has

two boundary layers_ The Green_s function is also unbounded as O______ at

x _ _ as _ _ __

__ Parabolic Problems

Bibliography

_ _ S_ Adjerid_ M_ Ai_a_ and J_E_ Flaherty_ Computational methods for singularly per_

turbed systems_ In J_D_ Cronin and Jr_ R_E_ O_Malley_ editors_ Analyzing Multiscale

Phenomena Using Singular Perturbation Methods_ volume __ of Proceedings of Sym_

posia in Applied Mathematics_ pages __%___ Providence_ ____ AMS_

___ U_M_ Ascher and L_R_ Petzold_ Computer Methods for Ordinary Dierential Equa_

tions and Dierential_Algebraic Equations_ SIAM_ Philadelphia_ ____

___ K_E_ Brenan_ S_L Campbell_ and L_R_ Petzold_ Numerical Solution of Initial_Value

Problems in Dierential_Algebraic Equations_ North Holland_ New York_ ____

___ J_E_ Flaherty_ A rational function approximation for the integration point in ex_

ponentially weighted _nite element methods_ International Journal of Numerical

Methods in Engineering_ _____%__ _ ____

___ J_E_ Flaherty and W_ Mathon_ Collocation with polynomial and tension splines for

singularly_perturbed boundary value problems_ SIAM Journal on Scie_nti_c and

Statistical Computation_ ____%____ ____

___ C_W_ Gear_ Numerical Initial Value Problems in Ordinary Dierential Equations_

Prentice Hall_ Englewood Cli_s_ __ _

___ E_ Hairer_ S_P_ Norsett_ and G_ Wanner_ Solving Ordinary Dierential Equations I_

Nonsti Problems_ Springer_Verlag_ Berlin_ second edition_ ____

___ E_ Hairer and G_ Wanner_ Solving Ordinary Dierential Equations II_ Sti and

Dierential Algebraic Problems_ Springer_Verlag_ Berlin_ __ _

___ C_ Johnson_ Numerical Solution of Partial Dierential Equations by the Finite Ele_

ment method_ Cambridge_ Cambridge_ ____

_ __ H__G_ Roos_ M_ Stynes_ and L_ Tobiska_ Numerical Methods for Singularly Perturbed

Dierential Equations_ Springer_Verlag_ Berlin_ ____

Chapter __

Пресс-релиз популярных книг

Parabolic Problems

Популярные книги

Популярные статьи

Навигация