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

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

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

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

Mesh Generation and Assembly

Back

___ Introduction

There are several reasons for the popularity of _nite element methods_ Large code seg_

ments can be implemented for a wide class of problems_ The software can handle complex

geometry_ Little or no software changes are needed when boundary conditions change_

domain shapes change_ or coe_cients vary_ A typical _nite element software framework

contains a preprocessing module to de_ne the problem geometry and data_ a processing

module to assemble and solve the _nite element system_ and a postprocessing module to

output the solution and calculate additional quantities of interest_ The preprocessing

module

_ creates a computer model of the problem domain __ perhaps_ using a computer

aided design _CAD system_

_ discretizes _ into a _nite element mesh_

_ creates geometric and mesh databases describing the mesh entities _vertices_ edges_

faces and elements and their relationships to each other and to the problem ge_

ometry_ and

_ de_nes problem_dependent data such as coe_cient functions_ loading_ initial data_

and boundary data_

The processing module

_ generates element stiness and mass matrices and load vectors_

_ assembles the global stiness and mass matrices and load vector_

_ enforces any essential boundary conditions_ and

_ Mesh Generation and Assembly

_ solves the linear _or nonlinear algebraic system for the _nite element solution_

The postprocessing modules

_ calculates additional quantities of interest_ such as stresses_ total energy_ and a

posteriori error estimates_ and

_ stores and displaying solution information_

In this chapter_ we study the preprocessing and processing steps with the exception of

the geometrical description and solution procedures_ The former topic is not addressed

while the latter subject will be covered in Chapter ___

___ Mesh Generation

Discretizing two_dimensional domains into triangular or quadrilateral _nite element meshes

can either be a simple or di_cult task depending on geometric or solution complexi_

ties_ Discretizing three_dimensional domains is currently not simple_ Uniform meshes

may be appropriate for some problems having simple geometric shapes_ but_ even there_

nonuniform meshes might provide better performance when solutions vary rapidly_ e_g__

in boundary layers_ Finite element techniques and software have always been associated

with unstructured and nonuniform meshes_ Early software left it to users to generate

meshes manually_ This required the entry of the coordinates of all element vertices_

Node and element indexing_ typically_ was also done manually_ This is a tedious and

error prone process that has largely been automated_ at least in two dimensions_ Adap_

tive solution_based mesh re_nement procedures concentrate meshes in regions of rapid

solution variation and attempt to automate the task of modifying _re_ningcoarsening

an existing mesh ___ __ __ __ ____ While we will not attempt a thorough treatment of

all approaches_ we will discuss the essential ideas of mesh generation by _i mapping

techniques where a complex domain is transformed into a simpler one where a mesh may

be easily generated and _ii direct techniques where a mesh is generated on the original

domain_

_____ Mesh Generation by Coordinate Mapping

Scientists and engineers have used coordinate mappings for some time to simplify ge_

ometric di_culties_ The mappings can either employ analytical functions or piecewise

polynomials as presented in Chapter __ The procedure begins with mappings

x _ f____ _ _ y _ f____ _

____ Mesh Generation _

that relate the problem domain in physical _x_ y space to its image in the simpler ___ _

space_ A simply connected region and its computational counterpart appear in Figure

______ It will be convenient to introduce the vectors

xT _ _x_ y__ f ___ _ T _ _f____ _ _ f____ _ _ ______a

and write the coordinate transformation as

x _ f ___ _ ______b

______

(

0,)

1,2

___

____

___

_________________________________________________________________

_____________

___________________________

____________

_ _

1,1

2,1

2,2

(

(



,0)

,1)

1,)





1,1 2,1

1,2 2,2

Figure ______ Mapping of a simply connected region _left onto a rectangular computa_

tional domain _right _

In Figure ______ we show a region with four segments f ___ _ _ f ___ _ _ f___ _ _ and f___ _

that are related to the computational lines _ _ __ _ _ __ _ _ __ and _ _ __ respectively_

_The four curved segments may involve dierent functions_ but we have written them all

as f for simplicity_

Also consider the projection operators

x _ P__f _ _N

___ f___ _ _ _N

___ f___ _ _ ______a

x _ P__f _ _N

___ f ___ _ _ _N

___ f ___ _ _ ______b

where

___ _ _ _ __ ______c

and

___ _ _ ______d

_ Mesh Generation and Assembly

are the familiar hat functions scaled to the interval _ _ _ _ __

As shown in Figure ______ the mapping x _ P__f transforms the left and right edges

of the domain correctly_ but ignores the top and bottom while the mapping x _ P__f

transforms the top and bottom boundaries correctly but not the sides_ Coordinate lines

of constant _ and _ are mapped as either curves or straight lines on the physical domain_

1,1

1,2

2,2 2,2

1,1

1,2

2,1 2,1

Figure ______ The transformations x _ P__f _left and x _ P__f _right as applied to

the simply_connected domain shown in Figure ______

1,1

2,1

2,2

1,2

1,1

2,1

1,2

2,2

Figure ______ Illustrations of the transformations x _ P_P__f _left and x _ P_ _P__f

_right as applied to the simply_connected domain shown in Figure ______

With a goal of constructing an eective mapping_ let us introduce the tensor product

and Boolean sums of the projections ______ as

x _ P_P__f _

Xi__

Xj__

i__ _N

j__ f _i _ __ j _ _ ______a

____ Mesh Generation _

x _ P_ _P__f _ P__f _ P__f _P_P__f _ ______b

An application of these transformations to a simply_connected domain is shown in Figure

______ The transformation ______a is a bilinear function of _ and _ while ______b is clearly

the one needed to map the simply connected domain onto the computational plane_ Lines

of constant _ and _ become curves in the physical domain _Figure _____ _

Although these transformations are simple_ they have been used to map relatively

complex two_ and three_dimensional regions_ Two examples involving the _ow about an

airfoil are shown in Figure ______ With the transformation shown at the top of the _gure_

the entire surface of the airfoil is mapped to _ _ _ ____ _ A cut is made from the trailing

edge of the airfoil and the curve so de_ned is mapped to the left __ _ __ ___ and right

__ _ __ ___ edges of the computational domain_ The entire far _eld is mapped to the top

__ _ __ ___ of the computational domain_ Lines of constant _ are rays from the airfoil

surface to the far _eld boundary in the physical plane_ Lines of constant _ are closed

curves encircling the airfoil_ Meshes constructed in this manner are called _O_grids__ In

the bottom of Figure ______ the surface of the airfoil is mapped to a portion ____ of the

_ axis_ The cut from the trailing edge is mapped to the rest ____ and ___ of the axis_

The _right out_ow boundary is mapped to the left ____ and right ____ edges of the

computational domain_ and the top_ left_ and bottom far _eld boundaries are mapped

to the top __ _ __ ___ of the computational domain_ Lines of constant _ become curves

beginning and ending at the out_ow boundary and surrounding the airfoil_ Lines of

constant _ are rays from the airfoil surface or the cut to the outer boundary_ This mesh

is called a _C_grid__

_____ Unstructured Mesh Generation

There are several approaches to unstructured mesh generation_ Early attempts used

manual techniques where point_coordinates were explicitly de_ned_ Semi_automatic mesh

generation required manual input of a coarse mesh which could be uniformly re_ned by

dividing each element edge into K segments and connecting segments on opposite sides of

an element to create K_ _triangular elements_ More automatic procedures use advancing

fronts_ point insertion_ and recursive bisection_ We_ll discuss the latter procedure and

brie_y mention the former_

With recursive bisection ____ a two_dimensional region _ is embedded in a square _uni_

verse_ that is recursively quartered to create a set of disjoint squares called quadrants_

Quadrants are related through a hierarchical quadtree structure_ The original square

universe is regarded as the root of the tree and smaller quadrants created by subdivi_

sion are regarded as ospring of larger ones_ Quadrants intersecting __ are recursively

_ Mesh Generation and Assembly

___________

______________________________________________________

_ _ _ _ _ _ _ _ _ _ _





2 3

1 4

airfoil

___________

_ _ _ _ _ _ _ _ _

_ _

_ _ _ _ _ _ _ _ _

_ _

_ _ _ _ _ _ _ _ _

_ _

_ _ _ _ _ _ _ _ _

_ _

_ _ _ _ _ _ _ _ _

_ _

_ _ _ _ _ _ _ _ _

_ _





airfoil

1 2 3 4

5 6

Figure ______ _O_grid_ _top and _C_grid_ _bottom mappings of the _ow about an airfoil_

quartered until a prescribed spatial resolution of _ is obtained_ At this stage_ quadrants

that are leaf nodes of the tree and intersect _ _ __ are further divided into small sets

of triangular or quadrilateral elements_ Severe mesh gradation is avoided by imposing a

maximal one_level dierence between quadrants sharing a common edge_ This implies a

maximal two_level dierence between quadrants sharing a common vertex_

A simple example involving a domain consisting of a rectangle and a region within a

curved arc_ as shown in Figure ______ will illustrate the quadtree process_ In the upper

portion of the _gure_ the square universe containing the problem domain is quartered

creating the one_level tree structure shown at the upper right_ The quadrant containing

the curved arc is quartered and the resulting quadrant that intersects the arc is quartered

again to create the three_level tree shown in the lower right portion of the _gure_ A

triangular mesh generated for this tree structure is also shown_ The triangular elements

are associated with quadrants of the tree structure_ Quadrants and a mixed triangular_

and quadrilateral_element mesh for a more complex example are shown in Figure ______

Elements produced by the quadtree and octree techniques may have poor geometric

shapes near boundaries_ A _nal _smoothing_ of the mesh improves element shapes and

____ Mesh Generation _

_ __ _

_ _ _

Boundary quadrant

Interior quadrant

Exterior quadrant

Finite element

Figure ______ Finite quadtree mesh generation for a domain consisting of a rectangle and

a region within a curved arc_ One_level _top and three_level _bottom tree structures

are shown_ The mesh of triangular elements associated with the three_level quadtree is

shown superimposed_

further reduces mesh gradation near ___ Element vertices on __ are moved along the

boundary to provide a better approximation to it_ Pairs of boundary vertices that are too

close to each other may be collapsed to a single vertex_ Interior vertices are smoothed by a

Laplacian operation that places each vertex at the _centroid_ of its neighboring vertices_

To be speci_c_ let i be the index of a node to be re_positioned_ xi be its coordinates_ Pi

be the set of indices of all vertices that are connected to Node i by an element edge_ and

Qi contain the indices of vertices that are in the same quadrant as Node i but are not

_ Mesh Generation and Assembly

Figure ______ Quadtree structure and mixed triangular_ and quadrilateral_element mesh

generated from it_

connected to it by an edge_ Then

xi _

_Pj_Pi

xj _Pj_Qi

_dim_Pi _ dim_Qi

______

where dim_S is the number of element vertices in set S_ Additional details appear in

____ Data Structures _

Baehmann et al_ ____

Arbitrarily complex two_ and three_dimensional domains may be discretized by quadtree

and octree decomposition to produce unstructured grids_ Further solution_based mesh

re_nement may be done by subdividing appropriate terminal quadrants or octants and

generating a new mesh locally_ This unites mesh generation and adaptive mesh re_ne_

ment by a common tree data structure ____ The underlying tree structure is also suitable

for load balancing on a parallel computer ___ ___

The advancing front technique constructs a mesh by _notching_ elements from __

and propagating this process into the interior of the domain_ An example is shown in

Figure ______ This procedure provides better shape control than quadtree or octree but

problems arise as the advancing fronts intersect_ L ohner ____ has a description of this and

other mesh generation techniques_ Carey ___ presents a more recent treatment of mesh

generation_

Figure ______ Mesh generation by the advancing front technique_

___ Data Structures

Unstructured mesh computation requires a data structure to store the geometric infor_

mation_ There is some ambiguity concerning the information that should be computed

at the preprocessing stage_ but_ at the very least_ the processing module would have to

know

_ the vertices belonging to each element_

_ the spatial coordinates of each vertex_ and

_ the element edges_ faces_ or vertices that are on ___

The processing module would need more information when adaptivity is performed_ It_

for example_ would need a link to the geometric information in order to re_ne elements

__ Mesh Generation and Assembly

along a curved boundary_ Even without adaptivity_ the processing software may want

access to geometric information when using elements with curved edges or faces _cf_

Section ___ _ If the _nite element basiswere known at the preprocessing stage_ space could

be reserved for edge and interior nodes or for a symbolic factorization of the resulting

algebraic system _cf_ Chapter __ _

Beall and Shephard ___ introduced a database and data structure that have great

_exibility_ It is suitable for use with high_order and hierarchical bases_ adaptive mesh

re_nement andor order variation_ and arbitrarily complex domains_ It has a hierarchical

structure with three_dimensional elements _regions having pointers to their bounding

faces_ faces having pointers to their bounding edges_ and edges having pointers to their

bounding vertices_ Those mesh entities _elements_ faces_ edges_ and vertices on domain

boundaries have pointers to relevant geometric structures de_ning the problem domain_

This structure_ called the SCOREC mesh database_ is shown in Figure ______ Nodes may

be introduced as _xed points in space to be associated with shape functions_ When done_

these may be located by pointers from any mesh entity_

Element

Face

Edge

Vertex

Geometric

Model

Entities

Figure ______ SCOREC hierarchical mesh database_

Let us illustrate the data structure for the two_dimensional domain shown in Figure

______ As shown in Figure ______ this mesh has __ faces _two_dimensional elements _ __

edges_ and __ vertices_ The face and edge_pointer information is shown in Table ______

Each edge has two pointers back to the faces that contain it_ These are shown within

brackets in the table_ The use of tables and integer indices for pointers is done for

convenience and does not imply an array implementation of pointer data_ The edge and

____ Data Structures __

vertex_pointer information and the vertex_point coordinate data are shown in Table ______

Backward pointers from vertices to edges and pointers from vertices and edges on the

boundary to the geometric database have not been shown to simplify the presentation_

We have shown a small portion of the pointer structure near Edge __ in Figure ______

Links between common entities allow the mesh to be traversed by faces_ edges_ or vertices

in two dimensions_ Problem and solution data is stored with the appropriate entities_

9 15

7 19

13 14

13 16

Figure ______ Example illustrating the SCOREC mesh database_ Faces are indexed

as shown at the upper left_ edge numbering is shown at the upper right_ and vertex

numbering is shown at the bottom_

__ Mesh Generation and Assembly

Face Edge Edge Edge

_ _ _ _ _ _ _ _ __ _ _ _ __

_ _ _ _ __ _ _ _ __ __ _ _ __

_ _ _ _ _ __ _ _ __ _ _ _ __

_ __ _ _ __ __ _ _ __ __ _ _ __

_ __ _ _ __ __ _ _ __ __ _ _ ___

_ _ _ _ _ __ _ _ __ __ _ _ __

_ __ _ _ __ _ _ _ _ __ _ _ ___

_ _ _ _ __ __ _ _ ___ __ _ _ _

__ __ ___ __ __ ___ ___ __ ___ ___

__ __ ___ ___ __ ___ ___ __ ___ ___

__ __ ___ _ _ __ ___ ___ __ ___ ___

__ __ ___ ___ __ ___ ___ __ ___ ___

__ __ ___ ___ __ ___ ___ __ ___ _

__ __ ___ _ __ ___ ___ __ ___ _

__ __ ___ _ __ ___ ___ __ ___ ___

__ __ ___ ___ __ ___ ___ __ ___ ___

__ __ ___ ___ __ ___ ___ __ ___ _

__ __ ___ ___ __ ___ _ __ ___ ___

Table ______ Face and edge_pointer data for the mesh shown in Figure ______ Backward

pointers from edges to their bounding faces are shown in brackets_

Face 14

Edge 18

Face 10 ...

...

To Edges 17 and 19 To Edges 24 and 23

To vertices 10 and 11

...

Figure ______ Pointer structure in the vicinity of Edge ___

The SCOREC mesh database contains more information than necessary for a typi_

cal _nite element solution_ For example_ the edge information may be eliminated and

faces may point directly to vertices_ This would be a more traditional _nite element

data structure_ Although it saves storage and simpli_es the data structure_ it may be

wise to keep the edge information_ Adaptive mesh re_nement procedures often work by

edge splitting and these are simpli_ed when edge data is available_ Edge information

also simpli_es the application of boundary conditions_ especially when the boundary is

____ Data Structures __

Edge Vertices Edge Vertices

_ _ _ __ _ __

_ _ _ __ __ __

_ _ _ __ __ _

_ _ _ __ __ __

__ _ _ __ __ __

__ _ __ __ __ __

__ __ __ __ _ __

Vertex Coordinates

_ _____ ____

_ ____ ____

_ _____ ____

_ ____ ____

_ _____ ____

__ _____ ____

__ ____ ____

__ ____ _____

__ ____ ____

Table ______ Edge and vertex_pointer data _left and vertex and coordinate data _right

for the mesh shown in Figure ______

curved_ Only pointers are required for the edge information and_ in many implementa_

tions_ pointers require less storage than integers_ Nevertheless_ let us illustrate face and

vertex information for the simple mesh shown in Figure ______ which contains a mixture

of triangular and quadrilateral elements_ The face_vertex information is shown in Table

_____ and the vertex_coordinate data is shown in Table ______ Assuming quadratic shape

functions on the triangles and biquadratic shape functions on the rectangles_ a traditional

data structure would typically add nodes at the centers of all edges and the centers of

the rectangular faces_ In this example_ the midside and face nodes are associated with

faces_ however_ they could also have been associated with vertices_

Without edge data_ the database generally requires additional a priori assumptions_

For example_ we could agree to list vertices in counterclockwise order_ Edge nodes could

follow in counterclockwise order beginning with the node that is closest in the coun_

terclockwise direction to the _rst vertex_ Finally_ interior nodes may be listed in any

order_ The choice of the _rst vertex is arbitrary_ This strategy is generally a compromise

between storing a great deal of data with fast access and having low storage costs but

having to recompute information_ We could further reduce storage_ for example_ by not

saving the coordinates of the edge nodes_

__ Mesh Generation and Assembly

1 2 3

4 6

9 10

11 13

16 19

21 22

(1) (2)

(3)

(4)

(5)

Figure ______ Sample _nite element mesh involving a mixture of quadratic approximations

on triangles and biquadratic approximations on rectangles_ Face indices are shown in

parentheses_

Face Vertices Nodes

_ _ _ _ _ _ __ __ __ __

_ _ _ _ _ __ __ __ __ __

_ _ _ _ __ __ __

Table ______ Simpli_ed face_vertex data for the mesh of Figure ______

The type of _nite element basis must also be stored_ In the present example_ we could

attach it to the face_vertex table_ With the larger database described earlier_ we could

attach it to the appropriate entity_ In the spirit of the shape function decomposition

described in Sections ___ and ____ we could store information about a face shape function

with the face and information about an edge shape function with the edge_ This would

allow us to use variable_order approximations _p_re_nement _

Without edge data_ we need a way of determining those edges that are on ___ This

can be done by adopting a convention that the edge between the _rst and second vertices

of each face is Edge __ Remaining edges are numbered in counterclockwise order_ A

sample boundary data table for the mesh of Figure _____ is shown on the right of Table

______ The _rst row of the table identi_es Edge _ of Face _ as being on a boundary of

the domain_ Similarly_ the second row of the table identi_es Edge _ of Face _ as being a

boundary edge_ etc_ Regions with curved edges would need pointers back to the geometric

database_

____ Coordinate Transformations __

Vertex Coordinates

_ ____ ____

__ ____ ____

Face Edge

_ _

Table ______ Vertex and coordinate data _left and boundary data _right for the _nite

element mesh shown in Figure ______

___ Coordinate Transformations

Coordinate transformations enable us to develop element stiness and mass matrices

and load vectors on canonical triangular_ square_ tetrahedral_ and cubic elements in a

computational domain and map these to actual elements in the physical domain_ Useful

transformations must _i be simple to evaluate_ _ii preserve continuity of the _nite

element solution and geometry_ and _iii be invertible_ The latter requirement ensures

that each point within the actual element corresponds to one and only one point in the

canonical element_ Focusing on two dimensions_ this requires the Jacobian

Je __ _ x_ x_

y_ y_ _ ______

of the transformation of Element e in the physical _x_ y _plane to the canonical element

in the computational ___ _ _plane to be nonsingular_

The most popular coordinate transformations are_ naturally_ piecewise_polynomial

__ Mesh Generation and Assembly

functions_ These mappings are called subparametric_ isoparametric_ and superparametric

when their polynomial degree is_ respectively_ lower than_ equal to_ and greater than

that used for the trial function_ As we have seen in Chapter __ the transformations use

the same shape functions as the _nite element solutions_ We illustrated linear _Section

___ and bilinear _Section ___ transformations for_ respectively_ mapping triangles and

quadrilaterals to canonical elements_ We have two tasks in front of us_ _i determining

whether higher_degree piecewise polynomial mappings can be used to advantage and _ii

ensuring that these transformations will be nonsingular_

Example ______ Recall the bilinear transformation of a _ _ _ canonical square to a

quadrilateral that was introduced in Section ___ _Figure _____

____________________

________

1,1 (x ,y )

(x ,y )

2,2 (x ,y )

1,2

(x ,y )







2,1

1,2 2,2

1,1 2,1

1 1



Figure ______ Bilinear mapping of a quadrilateral to a _ _ _ square_

_ x___ _

y___ _ _ _

Xi__

Xj__

_ xij

yij _Ni_j___ _ _ ______a

where

Ni_j___ _ _ _N

i__ _N

j__ _ i_ j _ __ __ ______b

and

i__ _ _ __ _ _ ___ if i _ _

__ _ _ ___ if i _ _

_ ______c

The vertices of the square ______ _ _____ _ ____ _ _ ___ _ are mapped to the vertices of

the quadrilateral _x____ y___ _ _x____ y___ _ _x____ y___ _ _x____ y___ _ The bilinear transforma_

tion is linear along each edge_ so the quadrilateral element has straight sides_

____ Coordinate Transformations __

Dierentiating ______a while using ______b_c

x_ _

x__ _ x__

___ _

x__ _ x__

___ _

y_ _

y__ _ y__

___ _

y__ _ y__

___ _

x_ _

x__ _ x__

___ _

x__ _ x__

___ _

y_ _

y__ _ y__

___ _

y__ _ y__

___ _

Substituting these formulas into ______ and evaluating the determinant reveals that the

quadratic terms cancel_ hence_ the determinant of Je is a linear function of _ and _ rather

than a bilinear function_ Therefore_ it su_ces to check that det_Je has the same sign at

each of the four vertices_ For example_

det_Je______ _ x_______ y_______ _ x_______ y_______

det_Je______ _ _x__ _ x__ _y__ _ y__ _ _x__ _ x__ _y__ _ y__ _

The cross product formula for two_component vectors indicates that

det_Je______ _ h_h_ sin ____

where h__ h__ and ___ are the lengths of two adjacent sides and the angle between them

_Figure _____ _ Similar formulas apply at the other vertices_ Therefore_ det_Je will not

vanish if and only if _ij _ at each vertex_ i_e__ if and only if the quadrilateral is convex_

Polynomial shape functions and bases are constructed on the canonical element as

described in Chapter __ For example_ the restriction of a bilinear _isoparametric trial

function to the canonical element would have the form

U___ _ _

Xi__

Xj__

ci_jNi_j___ _ _

A subparametric approximation might_ for example_ use a piecewise_bilinear coordinate

transformation ______ with a piecewise_biquadratic trial function_ Let us illustrate this

using the element node numbering of Section ___ as shown in Figure ______ Using ______ _

the restriction of the piecewise_biquadratic polynomial trial function to the canonical

element is

U___ _ _

Xi__

Xj__

ci_jN_

i_j___ _ ______a

__ Mesh Generation and Assembly

1,2

2,1

1,2 2,2

1,1 2,1

1,1

2,2

3,2

1,3

3,1

2,3

3,3 1,3 3,3 2,3

3,1

3,2





Figure ______ Bilinear mapping to a unit square with a biquadratic trial function_

_ _

__ _

2,1

1,2 2,2

1,1 2,1

1,1

2,2

1,3 3,3 2,3

3,1

3,2





1,3

3,2

3,1

2,3

3,3

1,2

Figure ______ Biquadratic mapping of the unit square to a curvilinear element_

where the superscript _ is used to identify biquadratic shape functions

i_j___ _ _ _N

i __ _N

j __ _ i_ j _ __ __ __ ______b

with

_N _

i __ _ __

____ _ _ ___ if i _ _

___ _ _ ___ if i _ _

_ _ ___ if i _ _

_ ______c

Example ______ A biquadratic transformation of the canonical square has the form

_ x___ _

y___ _ _ _

Xi__

Xj__

_ xij

yij _N_

i_j___ _ _ ______

where N_

i_j___ _ _ i_ j _ __ __ __ is given by ______ _

This transformation produces an element in the _x_ y _plane having curved _quadratic

edges as shown in Figure ______ An isoparametric approximation would be biquadratic

____ Coordinate Transformations __

4 2





Figure ______ Quadratic mapping of a triangle having one curved side_

and have the form of ______ _ The interior node ____ is awkward and can be eliminated

by using a second_order serendipity _cf_ Problems _____ or hierarchical transformation

_cf_ Section ___ _

Example ______ The biquadratic transformation described in Example _____ is useful

for discretizing domains having curved boundaries_ With a similar goal_ we describe a

transformation for creating triangular elements having one curved and two straight sides

_Figure _____ _ Let us approximate the curved boundary by a quadratic polynomial and

map the element onto a canonical right triangle by the quadratic transformation

_ x___ _

y___ _ _ _

Xi__

_ xi

yi _N_

i ___ _ _ ______a

where the quadratic Lagrange shape functions are _cf_ Problem _____

j _ _j_j _ ___ _ j _ __ __ __ ______b

_ _ ____ N_

_ _ ____ ______c

and

_ _ __ _ _ __ _ _ __ _ _ __ ______

Equations ______ and ______ describe a general quadratic transformation_ We have a

more restricted situation with

x_ _ _x_ _ x_ ___ y_ _ _y_ _ y_ ___

__ Mesh Generation and Assembly

This simpli_es the transformation ______a to

_ x___ _

y___ _ _ _ _ x_

y_ _ !N

_ _ _ x_

y_ _ !N

_ _ _ x_

y_ _ !N

_ _ _ x_

y_ _ !N

_ _ ______a

where_ upon use of ______ and ______ _

_ _ N_

_ _ _N_

_ _ N_

_ __ _ _ _ __ _ _ __ ______b

_ _ N_

_ __ _ ___ _ __ _ ______c

_ _ N_

_ __ _ ___ _ __ _ ______d

_ _ N_

_ _ ____ ______e

From these results_ we see that the mappings on edges ___ __ _ _ and ___ __ _ _ are

linear and are_ respectively_ given by

_ x

y _ _ _ x_

y_ ___ _ _ _ _ x_

y_ ___ _ x

y _ _ _ x_

y_ ___ _ _ _ _ x_

y_ ___

The Jacobian determinant of the transformation can vanish depending on the location

of Node __ The analysis may be simpli_ed by constructing the transformation in two

steps_ In the _rst step_ we use a linear transformation to map an arbitrary element onto

a canonical element having vertices at ___ _ _ ___ _ _ and ___ _ but with one curved side_

In the second step_ we remove the curved side using the quadratic transformation ______ _

The linear mapping of the _rst step has a constant Jacobian determinant and_ therefore_

cannot aect the invertibility of the system_ Thus_ it su_ces to consider the second step

of the transformation as shown in Figure ______ Setting _x__ y_ _ ___ _ _ _x__ y_ _ ___ _ _

and _x__ y_ _ ___ _ in ______a yields

_ x___ _

y___ _ _ _ _ _

_ _ !N

_ _ _ _

_ _ !N

_ _ _ x_

y_ _ !N

_ _

Using ______c_e

_ x___ _

y___ _ _ _ _ ___ _ __

___ _ __ __ ___ _ x_

y_ __

Calculating the Jacobian

Je___ _ _ _ x_ x_

y_ y_ _ _ _ _ _ __ _ _x__ ___ _ _x__

___ _ _y__ _ _ __ _ _y__ __

____ Coordinate Transformations __

________

4 2





1 4

Figure ______ Quadratic mapping of a right triangle having one curved side_ The shaded

region indicates where Node _ can be placed without introducing a singularity in the

mapping_

we _nd the determinant as

det_Je___ _ _ ____x_ _ _ _ _ __y_ _ _ __

The Jacobian determinant is a linear function of _ and __ thus_ as with Example ______

we need only ensure that it has the same sign at each of its three vertices_ We have

det_Je___ _ _ __ det_Je___ _ _ _x_ _ __ det_Je____ _ _y_ _ __

Hence_ the Jacobian determinant will not vanish and the mapping will be invertible when

x_ _ ___ and y_ _ ___ _cf_ Problem _ at the end of this section _ This region is shown

shaded on the triangle of Figure ______

Problems

__ Consider the second_order serendipity shape functions of Problem _____ or the

second_order hierarchical shape functions of Section ____ Let the four vertex nodes

be numbered ___ _ _ ___ _ _ ___ _ _ and ___ _ and the four midside nodes be numbered

___ _ _ ___ _ _ ___ _ _ and ___ _ _ Use the serendipity shape functions of Problem _____

to map the canonical _ _ _ square element onto an eight_noded quadrilateral ele_

ment with curved sides in the _x_ y _plane_ Assume that the vertex and midside

nodes of the physical element have the same numbering as the canonical element

but have coordinates at _xij_ yij _ i_ j _ __ __ __ i _ j __ __ Can the Jacobian of the

transformation vanish for some particular choices of _x_ y " _This is not a simple

question_ It su_ces to give some qualitative reasoning as to how and why the

Jacobian may or may not vanish_

__ Mesh Generation and Assembly

__ Consider the transformation ______ of Example _____ with x_ _ y_ _ ___ and sketch

the element in the _x_ y _plane_ Sketch the element for some choice of x_ _ y_ _ ____

___ Generation of Element Matrices and Vectors and

Their Assembly

Having discretized the domain_ the next step is to select a _nite element basis and

generate and assemble the element stiness and mass matrices and load vectors_ As a

review_ we summarize some of the two_dimensional shape functions that were developed

in Chapter _ in Tables _____ and ______ Nodes are shown on the mesh entities for the

Lagrangian and hierarchical shape functions_ As noted in Section ____ however_ the shape

functions may be associated with the entities without introducing modal points_ The

number of parameters np for an element having order p shape functions is presented for

p _ __ __ __ __ We also list an estimate of the number of unknowns _degrees of freedom N

for scalar problems solved on unit square domains using uniform meshes of _n_ triangular

or n_ square elements_

Both the Lagrange and hierarchical bases of order p have the same number of param_

eters and degrees of freedom on the uniform triangular meshes_ Without constraints for

Dirichlet data_ the number of degrees of freedom is N _ _pn __ _ _cf_ Problem _ at the

end of this section _ Dirichlet data on the entire boundary would reduce N by O_pn

and_ hence_ be a higher_order eect when n is large_ The asymptotic approximation

N _ _pn _ is recorded in Table ______ Similarly_ bi_polynomial approximations of order p

on squares with n_ uniform elements have N _ _pn _ _ _ degrees of freedom _again_ cf_

Problem _ _ The asymptotic approximation _pn _ is reported in Table ______ Under the

same conditions_ hierarchical bases on squares have

N _ _ __p _ _ n_ _ _pn _ __ if p _ _

_p_ _ p _ _ n___ _ _pn _ __ if p _

degrees of freedom_ The asymptotic values N _ __p _ _ N__ p _ __ and N _ _p_ _ p _

_ n____ p __ are reported in Table ______

The Lagrange and hierarchical bases on triangles and the Lagrange bi_polynomial

bases on squares have approximately the same number of degrees of freedom for a given

order p_ The hierarchical bases on squares have about half the degrees of freedom of the

others_ The bi_polynomial Lagrange shape functions on a square have the largest number

of parameters per element for a given p_ The number of parameters per element aects the

element matrix and vector computations while the number of degrees of freedom aects

the solution time_ We cannot_ however_ draw _rm conclusions about the superiority of

one basis relative to another_ The selection of an optimal basis for an intended level

____ Element Matrices and Their Assembly __

p Lagrange Hierarchical np N _ p_n_

Stencil Stencil

__ _ n_

_ _n_

__ _n_

__ __

__ __n_

Table ______ Shape function placement for Lagrange and hierarchical _nite element ap_

proximations of degrees p _ __ __ __ _ on triangular elements with their number of param_

eters per element np and degrees of freedom N on a square with _n_ elements_ Circles

indicate additional shape functions located on a mesh entity_

of accuracy is a complex issue that depends on solution smoothness_ geometry_ and

the partial dierential system_ We_ll examine this topic in a later chapter_ At least it

seems clear that bi_polynomial bases are not competitive with hierarchical ones on square

elements_

_____ Generation of Element Matrices and Vectors

The generation of the element stiness and mass matrices and load vectors is largely

independent of the partial dierential system being solved_ however_ let us focus on the

model problem of Section ___ in order to illustrate the procedures less abstractly_ Thus_

consider the two_dimensional Galerkin problem_ determine u H_

E satisfying

A_v_ u _ _v_ f _ _v H_

_ _ ______a

__ Mesh Generation and Assembly

p Lagrange Hierarchical

Stencil np N _ Mn_ Stencil np N _ Mn_

_ _ _

_ __ _

__ _ n_

_ _ _

_ __ _

__ _ n_

_ _ _

_ __ _

_ _ _

_ _n_

_ _ _

_ __ _

_ _ _

_ __ _

_ _n_

_ _ _

_ _

_ __ _

_ _ _

_ _

_ __ _

_ _ _

__ __ _n_

__ _n_

_ _ _

_ __ _

_ _ _

_ __ _

_ _ _

_ __ _

_ _ _

_ _ _ _ __ _

_ _ _

_ __ _

__ __n_

_ _ _

_ __ _

_ _ _

_ __ _

__ _n_

Table ______ Shape function placement for bi_polynomial Lagrange and hierarchical ap_

proximations of degrees p _ __ __ __ _ on square elements with their number of parameters

per element np and degrees of freedom N on a square with n_ elements_ Circles indicate

additional shape functions located on a mesh entity_

where

_v_ f _ ZZ

vfdxdy ______b

A_v_ u _ ZZ

_p_vxux _ vyuy _ qvu_dxdy ______c

As usual_ _ is a two_dimensional domain with boundary __ _ __E _ __N_ Recall that

smooth solutions of ______ satisfy

__pux x _ _puy y _ qu _ f_ _x_ y __ ______a

____ Element Matrices and Their Assembly __

u _ __ _x_ y __E_ ______b

un _ __ _x_ y __N_ ______c

where n is the unit outward normal vector to ___ Trivial natural boundary conditions

are considered for simplicity_ More complicated situations will be examined later in this

section_

Following the one_dimensional examples of Chapters _ and __ we select _nite_dimensional

subspaces SN

E and SN

_ of H_

E and H_

_ and write ______b_c as the sum of contributions

over elements

_V_ f _

Xe__

_V_ f e_ ______a

A_V_U _

Xe__

Ae_V_U _ ______b

Here_ N is the number of elements in the mesh_

_V_ f e _ ZZ

V fdxdy ______c

is the local L_ inner product_

Ae_V_U _ ZZ

_p_VxUx _ VyUy _ qV U_dxdy ______d

is the local strain energy_ and _e is the portion of _ occupied by element e_

The evaluation of ______c_d can be simple or complex depending on the functions

p_ q_ and f and the mesh used to discretize __ If p and q were constant_ for example_

the local strain energy ______d could be integrated exactly as illustrated in Chapters _

and _ for one_dimensional problems_ Let_s pursue a more general approach and discuss

procedures based on transforming integrals ______c_d on element e to a canonical element

_ and evaluating them numerically_ Thus_ let U____ _ _ U_x___ _ _ y___ _ and V____ _ _

V _x___ _ _ y___ _ and transform the integrals ______c_d to element _ to get

_V_ f e _ ZZ

V____ _ f_x___ _ _ y___ _ det_Je d_d__ ______a

Ae_V_U _ ZZ

_p_V__ _x _ V___x _U___x _ U___x _

__ Mesh Generation and Assembly

p_V___y _ V___y _U___y _ U___y _ qV_U__ det_Je d_d_

where Je is the Jacobian of the transformation _cf_ ______ _

Expanding the terms in the strain energy

Ae_V_U _ ZZ

_g_eV__U__ _ g_e_V__U__ _ V__U__ _ g_eV__U__ _ qV_U__ det_Je d_d_

______b

where

g_e _ p_x___ _ _ y___ _ ___

x _ __

y__ ______c

g_e _ p_x___ _ _ y___ _ __x_x _ _y_y__ ______d

g_e _ p_x___ _ _ y___ _ ___

x _ __

y__ ______e

The integrand of ______b might appear to be polynomial for constant p and a poly_

nomial mapping_ however_ this is not the case_ In Section ____ we showed that the inverse

coordinate mapping satis_es

_x _

det_Je

_ _y _ _

det_Je

_ _x _ _

det_Je

_ _y _

det_Je

_ ______

The functions gie_ i _ __ __ __ are proportional to ___ det_Je ___ thus_ the integrand of

______b is a rational function unless_ of course_ det_Je is a constant_

Let us write U_ and V_ in the form

U____ _ _ cTe

N___ _ _ N___ _ Tce_ V____ _ _ dTe

N___ _ _ N___ _ Tde ______

where the vectors ce and de contain the elemental parameters and N___ _ is a vector

containing the elemental shape functions_

Example ______ For a linear polynomial on the canonical right ___ triangular element

having vertices numbered _ to _ as shown in Figure ______

ce __

ce__

_ N___ _ __

_ _ _ _ _

The actual vertex indices_ shown as i_ j_ abd k_ are mapped to the canonical indices __

__ and __

Example ______ The treatment of hierarchical polynomials is more involved because

there can be more than one parameter per node_ Consider the case of a cubic hierarchical

____ Element Matrices and Their Assembly __

_ i

i, 1 j, 2

k, 3





Figure ______ Linear transformation of a triangular element e to a canonical right ___

triangle_

function on a triangle_ Translating the basis construction of Section ___ to the canonical

element_ we obtain an approximation of the form ______ with

cTe

_ _ce___ ce___ ____ ce____

NT _ _N____ _ _N____ _ _ ____N_____ _ __

The basis has ten shape functions per element _cf_ ________ _ which are ordered as

N____ _ _ _ _ __ _ _ __ N____ _ _ _ _ __ N____ _ _ _ _ __

N____ _ _ _

p____ N____ _ _ _

p____

N____ _ _ _

p_______ _ _ _ N____ _ _ _

p_______ _ _ _

N___ _ _ _

p______ _ __ _ N_____ _ _ ____

With this ordering_ the _rst three shape functions are associated with the vertices_ the

next three are quadratic corrections at the midsides_ the next three are cubic corrections

at the midsides_ and the last is a cubic _bubble function_ associated with the centroid

_Figure _____ _

An array implementation_ as described by ______ and Examples ____ _ and ______

may be the simplest data structure_ however_ implementations with structures linked to

geometric entities _Section ___ are also possible_

Substituting the polynomial representation ______ into the transformed strain energy

expression ______b and external load ______a yields

Ae_V_U _ dTe

_Ke _Me ce_ ______a

__ Mesh Generation and Assembly





1 4, 7 2

6, 9 5, 8

Figure ______ Shape function placement and numbering for a hierarchical cubic approxi_

mation on a canonical right ___ triangle_

_V_ f e _ dTe

fe_ ______b

where

Ke _ ZZ

_g_eN_NT_

_ g_e_N_NT_

_ N_NT_

_ g_eN_NT_

_ det_Je d_d__ ______a

Me _ ZZ

qNNT det_Je d_d__ ______b

fe _ ZZ

Nf det_Je d_d__ ______c

Here_ Ke and Me are the element stiness and mass matrices and fe is the element load

vector_ Numerical integration will generally be necessary to evaluate these arrays when

the coordinate transformation is not linear and we will study procedures to do this in

Chapter __

Element mass and stiness matrices and load vectors are generated for all elements

in the mesh and assembled into their proper locations in the global stiness and mass

matrix and load vector_ The positions of the elemental matrices and vectors in their

global counterparts are determined by their indexing_ In order to illustrate this point_

consider a linear shape function on an element with Vertices __ __ and _ as shown in

Figure ______ These vertex indices are mapped onto local indices_ e_g__ __ __ __ of the

canonical element and the correspondence is recorded as shown in Figure ______ After

generating the element matrices and vectors_ the global indexing determines where to add

____ Element Matrices and Their Assembly __

these entries into the global stines and mass matrix and load vector_ In the example

shown in Figure ______ the entry ke

__ is added to Row _ and Column _ of the global

stiness matrix K_ The entry ke

__ is added to Row _ and Column _ of K_ etc_

The assembly process avoids the explicit summations implied by ______ and yields

A_V_U _ dT _K _M c_ ______a

_V_ f _ dT f _ ______b

where

cT _ _c__ c__ ____ cN__ ______c

dT _ _d__ d__ ____ dN__ ______d

where K is the global stiness matrix_ M is the global mass matrix_ f is the global load

vector_ and N is the dimension of the trial space _or the number of degrees of freedom _

Imposing the Galerkin condition ______a

A_V_U _ _V_ f _ dT __K _M c _ f_ _ __ _d _N_ _______a

yields

_K_M c _ f _ _______b

_____ Essential and Neumann Boundary Conditions

It_s customary to ignore any essential boundary conditions during the assembly phase_

Were boundary conditions not imposed_ the matrix K_M would be singular_ Essential

boundary conditions constrain some of the ci_ i _ __ __ ____N_ and they must be imposed

before the algebraic system _______b can be solved_ In order to simplify the discussion_

let us suppose that either M _ _ or that M has been added to K so that _______ may

be written as

dT _Kc _ f_ _ __ _d _N_ _______a

Kc _ f _ _______b

__ Mesh Generation and Assembly

4 7

Global Local

_ _

1 2

Ke __

__ ke

fe __

_ _ _ _ _ _ _ _ _

K _

______________

_ke

__ _ke

_ke

__ _ke

_ke

__ _ke

____________

f _

____________

_fe

__________

Figure ______ Assembly of an element stiness matrix and load vector into their global

counterparts for a piecewise_linear polynomial approximation_ The actual vertex indices

are recorded and stored _top _ the element stiness matrix and load vector are calculated

_center _ and the indices are used to determine where to add the entries of the elemental

matrix and vector into the global stiness and mass matrix_

____ Element Matrices and Their Assembly __

Essential boundary conditions may either constrain a single ci or impose constraints

between several nodal variables_ In the former case_ we partition _______a as

_d_d__ __ K__ K__

K__ K__ __ c_

c_ _ _ _ f_

f_ _ _ __ _______a

where the essential boundary conditions are

c_ _ ___ _______b

Recall _Chapters _ and _ _ that the test function V should vanish on __E_ thus_ corre_

sponding to _______b

d_ _ __ _______c

The second _block_ of equations in _______a should never have been generated and_

actually_ we should have been solving

_ _K__c_ _ K__c_ _ f__ _ dT

_ _K__c_ _ K____ _ f__ _ __ _______a

Imposing the Galerkin condition that _______a vanish for all d__

K__c_ _ f_ _ K_____ _______b

Partitioning _______ need not be done explicitly as in _______ _ It can be done im_

plicitly without rearranging equations_ Consider the original system _______b

_______

k__ k_j k_N

___

kj_ kjj kjN

___

kN_ kNj kNN

_____

_______

___

_____

_______

___

_____

_ _______

Suppose that one boundary condition speci_es cj _ _j_ then the j th equation _row of

the system is deleted_ cj is replaced by the boundary condition_ and the coe_cients of cj

are moved to the right_hand side to obtain

_________

k__ k__j__ k__j__ k_N

___

kj____ kj___j__ kj___j__ kj___N

___

kN_ kN_j__ kN_j__ kNN

_______

_________

___

cj__

___

_______

_________

f_ _ k__j_j

___

fj__ _ kj___j_j

___

fN _ kN_j_j

_______

__ Mesh Generation and Assembly

When the algebraic system is large_ the cost of moving data when rows and columns

are removed from the system may outweigh the cost of solving a larger algebraic system_

In this case_ the boundary condition cj _ _j can be inserted as the j th equation of

_______ _ Although not necessary_ the j th column is usually moved to the right_hand

side to preserve symmetry_ The resulting larger problem is

_______

k__ _ k_N

___

_ _ _

___

kN_ _ kNN

_____

_______

___

_____

_______

f_ _ k__j_j

___

fN _ kN_j_j

_____

The treatment of essential boundary conditions that impose constraints among several

nodal variables is much more di_cult_ Suppose_ for example_ there are l boundary

conditions of the form

Tc _ __ _______

where T is an l_N matrix and _ is an l_vector_ In vector systems of partial dierential

equations_ such boundary conditions arise when constraints are speci_ed between dierent

components of the solution vector_ In scalar problems_ conditions having the form _______

arise when a _global_ boundary condition like

uds _ _

is speci_ed_ They could also arise with periodic boundary conditions which might_ for

example_ specify u___ y _ u___ y if u were periodic in x on a rectangle of unit length_

One could possibly solve _______ for l values of ci_ i _ __ __ ____N_ in terms of the

others_ Sometimes there is an obvious choice_ however_ often there is no clear way to

choose the unknowns to eliminate_ A poor choice can lead to ill_conditioning of the

algebraic system_ An alternate way of treating problems with boundary conditions such

as _______ is to embed Problem _______ in a constrained minimization problem which

may be solved using Lagrange multipliers_ Assuming K to be symmetric and positive

semi_de_nite_ _______ can be regarded as the minimum of

I_c_ _ cTKc _ _cT f _

Using Lagrange multipliers_ we minimize the modi_ed functional

#I_c_ __ _ cTKc _ _cT f _ __T _Tc _ _ _

____ Element Matrices and Their Assembly __

where _ is an l_vector of Lagrange multipliers_ Minimizing #I with respect to c and _

yields

_ K TT

T _ __ c

_ _ _ _ f

_ __ _______

The system _______ may or may not be simple to solve_ If K is non_singular then the

algorithm described in Problem _ at the end of this section is eective_ However_ since

boundary conditions are prescribed by _______ _ K may not be invertible_

Nontrivial Neumann boundary conditions on __N require the evaluation of an extra

line integral for those elements having edges on __N_ Suppose_ for example_ that the

variational principle ______ is replaced by_ determine u H_

E satisfying

A_v_ u _ _v_ f _ _ v__ __ _v H_

_ _ _______a

where

_ v__ __ Z_ N

v__x_ y ds_ _______b

s being a coordinate on __N_ As discussed in Chapter __ smooth solutions of _______

satisfy ______a _ the essential boundary conditions ______b _ and the natural boundary

condition

pun _ __ _x_ y __N_ _______

The line integral _______b is evaluated in the same manner as the area integrals and it

will alter the load vector f _cf_ Problem _ at the end of this section _

Problems

__ Determine the number of degrees of freedom when a scalar _nite element Galerkin

problem is solved using either Lagrange or hierarchical bases on a square region

having a uniform mesh of either _n_ triangular or n_ square elements_ Express

your answer in terms of p and n and compare it with the results of Tables _____

and ______

__ Assume that K is invertible and show that the following algorithm provides a

solution of _______ _

Solve KW _ TT for W

Let Y _ TW

Solve Ky _ f for y

Solve Y_ _ Ty_ _ for _

Solve Kc _ f _ TT_ for c

__ Mesh Generation and Assembly

__ Calculate the eect on the element load vector fe of a nontrivial Neumann condition

having the form _______ _

__ Consider the solution of Laplace_s equation

uxx _ uyy _ __ _x_ y __

on the unit square _ __ f_x_ y j_ _ x_y _ _g with Dirichlet boundary conditions

u _ __ _x_ y ___

As described in the beginning of this section_ create a mesh by dividing the unit

square into n_ uniform square elements and then into _n_ triangles by cutting each

square element in half along its positive sloping diagonal_

____ Using a Galerkin formulation with a piecewise_linear basis_ develop the element

stiness matrices for each of the two types of elements in the mesh_

____ Assemble the element stiness matrices to form the global stiness matrix_

____ Apply the Dirichlet boundary conditions and exhibit the _nal linear algebraic

system for the nodal unknowns_

__ The task is is to solve a Dirichlet problem on a square using available _nite element

software_ The problem is

_uxx _ uyy _ f_x_ y _ __ _x_ y __

with u _ _ on the boundary of the unit square _ _ f_x_ y j_ _ x_ y _ _g_ Select

f_x_ y so that the exact solution of the problem is

u_x_ y _ exy sin x sin _ y_

The Galerkin form of this problem is to _nd u H_

_ satisfying

_vxux _ vyuy _ vf_dxdy _ __ _v H_

_ _

Solve this problem on a sequence of _ner and _ner grids using piecewise linear_

quadratic_ and cubic _nite element bases_ Select a basic grid with either two or

four elements in it and obtain _ner grids by uniform re_nement of each element

into four elements_ Present plots of the energy error as functions of the number of

degrees of freedom _DOF _ the mesh spacing h_ and the CPU time for the three

polynomial bases_ De_ne h as the square root of the area of an average element_

____ Element Matrices and Their Assembly __

You may combine the convergence histories for the three polynomial solutions on

one graph_ Thus_ you_ll have three graphs_ error vs_ h_ error vs_ DOF_ and error

vs_ CPU time_ each having results for the three polynomial solutions_ Estimate the

convergence rates of the solutions_ Comment on the results_ Are they converging

at the theoretical rates" Are there any unexpected anomalies" If so_ try to explain

them_ You may include plots of solutions andor meshes to help answer these

questions_

__ Consider the Dirichlet problem for Laplace_s equation

$u _ uxx _ uyy _ __ _x_ y __

u_x_ y _ __x_ y _ _x_ y ___

where _ is the L_shaped region with lines connecting the Cartesian vertices ____ _

____ _ ____ _ _____ _ ______ _ _____ _ ____ _ Select __x_ y so that the exact solution

expressed in polar coordinates is

u_r_ _ r___ sin

with

x _ r cos _ y _ r sin _

This solution has a singularity in its _rst derivative at r _ __ The singularity

is typical of those associated with the solution of elliptic problems at re_entrant

corners such as the one found at the origin_

Because of symmetries_ the problem need only be solved on half of the L_shaped

domain_ i_e__ the trapezoidal region #_

with lines connecting the Cartesian vertices

____ _ ____ _ ____ _ _____ _ ____ _

The Galerkin form of this problem consists of determining u H_

_ _

vxux _

vyuy_

dxdy _

_ _

Functions u H_

E satisfy the essential boundary conditions

u_x_ y _ __ y _ __ _ _ x _ __

u_r_ _ r___ sin

_ x _ __ _ _ y _ __ y _ __ __ _ x _ __

These boundary conditions may be expressed in Cartesian coordinates by using

r_ _ x_ _ y__ tan _

__ Mesh Generation and Assembly

The solution of the Galerkin problem will also satisfy the natural boundary condi_

tion un _ u_ _ _ along the diagonal y _ _x_

Solve this problem using available _nite element software_ To begin_ create a

three_element initial mesh by placing lines between the vertices ____ and ____

and between ____ and ____ _ Generate _ner meshes by uniform re_nement and use

piecewise_polynomial bases of degrees one through three_

As in Problem __ present plots of the energy error as functions of the number of

degrees of freedom_ the mesh spacing h_ and the CPU time for the three polyno_

mial bases_ You may combine the convergence histories for the three polynomial

solutions on one graph_ De_ne h as the square root of the area of an average ele_

ment_ Estimate the convergence rates of the solutions_ Is accuracy higher with a

high_order method on a coarse mesh or with a low_order method on a _ne mesh"

If adaptivity is available_ use a piecewise_linear basis to calculate a solution using

adaptive h_re_nement_ Plot the energy error of this solution with those of the

uniform_mesh solutions_ Is the adaptive solution more e_cient" E_ciency may

be de_ned as less CPU time or fewer degrees of freedom for the same accuracy_

Contrast the uniform and adaptive meshes_

___ Assembly of Vector Systems

Vector systems of partial dierential equations may be treated in the same manner as the

scalar problems described in the previous section_ As an example_ consider the vector

version of the model problem ______ _ determine u H_

E satisfying

A_v_ u _ _v_ f _ _v H_

_ _ ______a

where

_v_ f _ ZZ

vT fdxdy_ ______b

A_v_ u _ ZZ

_vT

xPux _ vT

y Puy _ vTQu_dxdy_ ______c

The functions u_x_ y _ v_x_ y _ and f _x_ y are m_vectors and P and Q are m_m matrices_

Smooth solutions of ______ satisfy

__Pux x _ _Puy y _ Qu _ f _ _x_ y __ ______a

____ Assembly of Vector Systems __

u _ __ _x_ y __D_ un _ __ _x_ y __N_ ______b

Example ______ Consider the biharmonic equation

$_w _ f_x_ y _ _x_ y __

where

$_ __ _ xx _ _ yy

is the Laplacian and _ is a bounded two_dimensional region_ Problems involving bihar_

monic operators arise in elastic plate deformation_ slow viscous _ow_ combustion_ etc_

Depending on the boundary conditions_ this problem may be transformed to a system of

two second_order equations having the form ______ _ For example_ it seems natural to let

u_ _ _$w

then

_$u_ _ f_

Let w _ _u_ to obtain the vector system

_$u_ _ f_ _$u_ _ u_ _ __ _x_ y __

This system has the form ______ with

u _ _ u_

u_ __ P _ I_ Q _ _ _ _

_ _ __ f _ _ f

_ __

The simplest boundary conditions to prescribe are

w _ _u_ _ ___ $w _ _u_ _ ___ _x_ y ___

With these _Dirichlet boundary conditions_ the variational form of this problem is

______a with

_v_ f _ ZZ

v_fdxdy

and

A_v_ u _ ZZ

__v_ x_u_ x _ _v_ y_u_ y _ _v_ x_u_ x _ _v_ y_u_ y _ v_u__dxdy_

The requirement that ______a be satis_ed for all vector functions v H_

_ gives the two

scalar variational problems

__v_ x_u_ x _ _v_ y_u_ y _ v_f_dxdy _ __ _v_ H_

_ _

__ Mesh Generation and Assembly

__v_ x_u_ x _ _v_ y_u_ y _ v_u__dxdy _ __ _v_ H_

_ _

We may check that smooth solutions of these variational problems satisfy the pair of

second_order dierential equations listed above_

We note in passing_ that the boundary conditions presented with this example are

not the only ones of interest_ Other boundary conditions do not separate the system as

neatly_

Following the procedures described in Section ____ we evaluate ______ in an element_

by_element manner and transform the elemental strain energy and load vector to the

canonical element to obtain

Ae_V_U _ ZZ

_VT

G_eU__ _ VT

G_eU___

G_eU__ _ VT

_QU__ det_Je d_d__ ______a

where

G_e _ P___

x _ __

y __ G_e _ P__x_x _ _y_y__ G_e _ P___

x _ __

y __ ______b

and

_V_ f e _ ZZ

_ f det_Je d_d__ ______c

The restriction of the piecewise_polynomial approximation U_ to element e is written

in terms of shape functions as

U____ _ _

Xj__

ce_jNj___ _ ______a

where np is the number of shape functions on element e_ We have divided the vector

ce of parameters into its contributions ce_j _ j _ __ __ ____ n_ from the shape functions of

element e_ Thus_ we may write

cTe

_ _cT

e___ cT

e___ ____ cT

e_np__ ______b

In this form_ we may write U_ as

U_ _ NT ce _ cTe

N ______c

____ Assembly of Vector Systems __

where N is the npm _ m matrix

NT _ _N_I_N_I_ ____NnpI__ ______d

and the identity matrices have the dimension m of the partial dierential system_ The

simple linear shape functions will illustrate the formulation_

Example ______ Consider the solution of a system of m _ _ equations using a

piecewise_linear _nite element basis on triangles_ Suppose_ for convenience_ that the

node numbers of element e are __ __ and __ In order to simplify the notation_ we suppress

the subscript _ on U_ and V_ and the subscript e on ce_ The linear approximation on

element e then takes the form

_ U_

U_ _ _ _ c__

c__ _N____ _ _ _ c__

c__ _N____ _ _

where

N____ _ _ _ _ _ _ __ N____ _ _ __ N____ _ _ __

The _rst subscript on cij denotes its index in c and the second subscript identi_es the

vertex of element e_ The expression ______c takes the form

_ U_

U_ _ _ _ N_ _ N_ _ N_ _

_ N_ _ N_ _ N_ _

________

c__

______

Substituting ______c and a similar expression for V_ into ______a_e yields

Ae_V_U _ dTe

_Ke _Me ce_ _V_ f e _ dTe

fe_ ______

where

Ke _ ZZ

_N_G_eNT

_ N_G_eNT

_ det_Je d_d__ ______a

Me _ ZZ

NQNT det_Je d_d__ fe _ ZZ

Nf det_Je d_d__ ______b

Problems

__ Mesh Generation and Assembly

__ It is_ of course_ possible to use dierent shape functions for dierent solution com_

ponents_ This is often done with incompressible _ows where the pressure is approx_

imated by a basis having one degree less than that used for velocity_ Variational

formulations with dierent _elds are called mixed variational principles_ The result_

ing _nite element formulations are called mixed methods_ As an example_ consider

a vector problem having two components_ Suppose that a piecewise_linear basis is

used for the _rst variable and piecewise quadratics are used for the second_ Using

hierarchical bases_ select an ordering of unknowns and write the form of the _nite el_

ement solution on a canonical two_dimensional element_ What are the components

of the matrix N" For this approximation_ develop a formula for the element sti_

ness matrix ______a _ Express your answer in terms of the matrices Gie_ i _ __ __ __

and integrals of the shape functions_

Bibliography

___ I_ Babu%ska_ J_E_ Flaherty_ W_D_ Henshaw_ J_E_ Hopcroft_ J_E_ Oliger_ and T_ Tez_

duyar_ editors_ Modeling_ Mesh Generation_ and Adaptive Numerical Methods for

Partial Di erential Equations_ volume __ of The IMA Volumes in Mathematics and

its Applications_ New York_ _____ Springer_Verlag_

___ P_L_ Baehmann_ M_S_ Shephard_ and J_E_ Flaherty_ Adaptive analysis for automated

_nite element modeling_ In J_R_ Whiteman_ editor_ The Mathematics of Finite

Elements and Applications VI_ MAFELAP ____ pages ___&____ London_ _____

Academic Press_

___ P_L_ Baehmann_ S_L_ Witchen_ M_S_ Shephard_ K_R_ Grice_ and M_A_ Yerry_ Ro_

bust_ geometrically_based_ automatic two_dimensional mesh generation_ Interna

tional Journal of Numerical Methods in Engineering_ _______&_____ _____

___ M_W_ Beall and M_S_ Shephard_ A general topology_based mesh data structure_

International Journal of Numerical Methods in Engineering_ _______&_____ _____

___ M_W_ Bern_ J_E_ Flaherty_ and M_ Luskin_ editors_ Grid Generation and Adaptive

Algorithms_ volume ___ of The IMA Volumes in Mathematics and its Applications_

New York_ _____ Springer_

___ G_F_ Carey_ Computational Grids_ Generation_ Adaptation_ and Solution Strategies_

Series in Computational and Physical Processes in Mechanics and Thermal science_

Taylor and Francis_ New York_ _____

___ J_E_ Flaherty_ R_ Loy_ M_S_ Shephard_ B_K_ Szymanski_ J_ Teresco_ and L_ Ziantz_

Adaptive local re_nement with octree load_balancing for the parallel solution of

three_dimensional conservation laws_ Parallel and Distributed Computing_ _____ to

appear_

___ J_E_ Flaherty_ R_M_ Loy_ C_ Ozturan_ M_S_ Shephard_ B_K_ Szymanski_ J_D_ Teresco_

and L_H_ Ziantz_ Parallel structures and dynamic load balancing for adaptive _nite

element computation_ Applied Numerical Mathematics_ ______&____ _____

__ Mesh Generation and Assembly

___ J_E_ Flaherty_ P_J_ Paslow_ M_S_ Shephard_ and J_D_ Vasilakis_ editors_ Adaptive

methods for Partial Di erential Equations_ Philadelphia_ _____ SIAM_

____ R_ L ohner_ Finite element methods in CFD_ Grid generation_ adaptivity and par_

allelization_ In H_ Deconinck and T_ Barth_ editors_ Unstructured Grid Methods for

Advection Dominated Flows_ number AGARD Report AGARD_R_____ Neuilly sur

Seine_ _____ Chapter __

____ R_ Verf urth_ A Review of Posteriori Error Estimation and Adaptive Mesh

Re_nement Techniques_ Teubner_Wiley_ Stuttgart_ _____

Chapter _

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

Mesh Generation and Assembly

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

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

Навигация