Abstract:
It provides a finite element method library which improves the reliability of a program using the finite element method as a library, and avoids calculation errors and an increase in convergence time due to programming errors. To this end, a library that describes a program process based on the finite element method is characterized in that a vector of a vector space spanned by basis functions of the finite element method, and a dual vector of a dual vector space defined by a metric derived from an inner product which is determined by the square integrations of the basis functions, are defined as different abstract data types.

Description:
FIELD OF THE INVENTION  
         [0001]    The present invention relates to the configuration of a basic simulator library of a finite element method.  
         BACKGROUND OF THE INVENTION  
         [0002]    According to Iwanami&#39;s Dictionary of Mathematics (third edition), discovery of the finite element method is accredited to either the article of Courant (Bull. Amer. Math. Soc., 49 (1943)) or the article of M. J. Turner et al., (J. Aero. Sci. 23 (1956) 805-823). In either case, the contribution of that method to the engineering field is still enormous almost half a century after its discovery; many articles are published every year and many textbooks are written. For example, O. C. Zienkiewicz, Masao Yoshiki &amp; Yoshiaki Yamada (TR), “Matrix Finite Element Method”, third revised edition, (Baifukan, 1984), and Motoki Yagawa &amp; Shinobu Yoshimura, “Finite Element Method” (Baifukan, 1991) are known.  
           [0003]    Furthermore, inventions that pertain to the finite element methods have been made, and Japanese Patent Publication No. 61-10771 is known as an early one. Especially, various proposals have been made in association with a mesh generation method. For example, Japanese Patent Laid-Open No. 6-231216, 6-231217, 9-185729, 9-259300, and the like are known.  
           [0004]    The basic concept of the finite element method grounds on the Ritz&#39;s method that obtains the solution of a partial differential equation defined as a boundary-value problem by approximating a function space F as a domain of definition of differential operators by a real finite dimensional Hilbert space H (since a real matrix expression is available even when the function space F assumes complex number values, only a case of a real number is explained without losing generality). As the inner product of the Hilbert space H, a natural inner product derived from the square integral (L2) norm (Iwanami&#39;s Dictionary of Mathematics, third edition, Iwanami, p. 61, 1994) is used. Also, as in the Ritz&#39;s method, approximation is evaluated in a weak topology. The finite dimensional Hilbert space H is defined by a set (V, &lt;,&gt;) of a finite dimensional vector space V and inner products &lt;,&gt;. In this case, linear mapping from V to F is designed to be injection, and various operations (e.g., differential operations) to F can be naturally limited to V as those operations to V.  
           [0005]    Many finite element methods are characterized by using a function which has a very narrow support (closure of a nonzero set) as a basis function of the finite dimensional Hilbert space. For this purpose, a metric M determined by natural inner products &lt;,&gt;; in the Hilbert space (V, &lt;,&gt;) with respect to the finite dimensional vector space V, i.e., linear mapping from V to its dual space V*  
           M:V→V*  (1)  
           [0006]    is generally obtained as a sparse matrix. Note that M is determined so that an inner product  
           &lt;,&gt;:V×V→R  (2)  
           [0007]    as bi-linear mapping to real number R matches  
           (,M):V×V→R, ((u,Mv)ε E R),  (3)  
           [0008]    using natural pairing (,):V×V*→R. Also, differential operators to the function space F can be considered as bounded operators (i.e., having real number coefficient finite dimensional matrix listing) on the finite dimensional Hilbert space H by the aforementioned natural limitations. Since the above limitations are not always homomorphic as an operator ring, they are not a natural inclusion as an operator ring. However, it is important to consider the vector space V as a natural inclusion of the function space F as a real number value vector space, and it is also important to consider that the finite element method is a numerical calculus technique defined on a finite dimensional operator ring.  
           [0009]    Since a function having a very narrow support is used as a basis function of V, as described above, a matrix expression of metric M and many physical differential operators in H normally becomes a very sparse matrix (a matrix which has a very small ratio of the number of nonzero elements to the total degree of freedom). For this reason, numerical calculations can be very easily made.  
           [0010]    In this manner, the mechanism itself of the finite element method has universality and versatility. An automatic generation program for automatically generating a program of the finite element method, and a finite element method library that provides a basic part of the finite element method as a library have been developed using such universality and versatility.  
           [0011]    As examples of the automatic generation program, Japanese Patent Publication Nos. 7-120275 and 7-120276, and Patent No. 2765888 are known. However, the automatic generation program has little flexibility since it is not easy to change inputs/outputs in correspondence with subjects.  
           [0012]    On the other hand, as an example of the finite element method library, a library named Diffpack is introduced in Hans Petter Langtangen, “Computational Partial Differential Equations: Numerical Methods and Diffpack Programming, (Lecture Notes in Computational Science and Engineering, 2, Springer, 1999. Furthermore, formulation of the finite element method on the hierarchical structure lattices shown in FIG. 3 as a base for a substantial finite element method library has been proposed. For example, such proposal has been described in detail in Michael Griebel &amp; Gerhard Zumbusch, “Parallel adaptive subspace correction schemes with applications to elasticity”, Computer Methods in applied mechanics and engineering 184, p. 303-332 (2000) and its references.  
           [0013]    Note that the finite element method library includes fundamental variable declarations, and subroutine and function sets, which are used to complete a self-produced finite element method program corresponding to each individual problem by defining variables in accordance with each individual problem and combining functions and the like provided by the library by the programmer as the user. Since processes by the finite element methods are mathematically formulated in the library, the method using the library can easily generate a finite element method program according to each individual problem (e.g., the strength of materials, mechanics of continua, electromagnetics, and the like) compared to generation of a finite element method program using the aforementioned automatic generation program.  
           [0014]    However, the user often uses the finite element method as a library without recognizing the metric given by formula (1). Many textbooks of linear algebra normally are based on an orthogonal normal basis as a base of a vector space, and do not consider a dual space and metric. As a result, the user introduces an inner product which is different from a natural inner product of the finite element method in the vector space, and such use distorts the original Hilbert space H of the finite element method, resulting in a quite different Hilbert space. The inner product of the Hilbert space is (bi-linear) mapping from a vector space to real numbers as quantities (ordered set), the magnitudes of which are to be compared, and the change of this inner product means definition of a quite different magnitude relation. In practice, since it is natural to make error evaluation or the like using this mapping to real numbers, the change in inner product often results in serious errors.  
           [0015]    In general, an operator ring is defined as a ring of operators on the Hilbert space, and the inner product of the Hilbert space plays a very important role in the operator ring. In case of a function space as an infinite dimensional vector space, a target phenomenon is often largely converted depending on the inner product. As described above, it is natural to consider the finite element method as the finite dimensional operator ring. However, since the finite element method simulates a physical space as an infinite dimensional vector space, the inner product must be carefully taken account. Hence, from such point of view, careless use of the inner product results in loss of mathematical rationality.  
           [0016]    In fact, it is a common practice in the finite element method to define, as the basis function, a connected support (closure of a domain of definition having nonzero values) on a cell or on the vicinity of that cell of non-structured lattices shown in FIG. 4. Individual cells have different volumes. On the other hand, nonzero matrix elements of a metric matrix have strong correlation to the volumes of individual cells in the finite element method, and a large volume ratio of elements means a large difference among the matrix elements. Hence, definition of another inner product without respect to the metric M results in delay of convergence and deterioration of precision. This poses problems of precision and reliability of the finite element method library itself, and is an important subject upon generating a library.  
           [0017]    In this manner, it is important to generate a program that can correctly calculate the inner product in the Hilbert space upon using the finite element method as a library in terms of the reliability of a program such as a reduction of calculation errors, an improvement in convergence time, and the like.  
         SUMMARY OF THE INVENTION  
         [0018]    The present invention has been made in consideration of the above problems, and has as its object to provide a finite element method library which can improve the reliability of a program using the finite element method as a library, and can avoid calculation errors and deterioration of a convergence time due to programming errors.  
           [0019]    In order to achieve the above object, for example, a finite element method library of the present invention comprises the following arrangement. That is, there is provided a finite element method library that describes a program process based on a finite element method, wherein a vector of a vector space spanned by basis functions of the finite element method, and a dual vector of a dual vector space defined by a metric derived from an inner product which is determined by square integrations of the basis functions, are defined as different abstract data types.  
           [0020]    Other features and advantages of the present invention will be apparent from the following description taken in conjunction with the accompanying drawings, in which like reference characters designate the same or similar parts throughout the figures thereof.  
       
    
    
     BRIEF DESCRIPTION OF THE DRAWINGS  
       [0021]    The accompanying drawings, which are incorporated in and constitute a part of the specification, illustrate embodiments of the invention and, together with the description, serve to explain the principles of the invention.  
         [0022]    [0022]FIG. 1 is a block diagram showing an example of a processing apparatus which comprises a finite element method library according to an embodiment of the present invention;  
         [0023]    [0023]FIG. 2 is a view showing the sequence in which a programmer programs a finite element method program using the processing apparatus according to the embodiment of the present invention;  
         [0024]    [0024]FIG. 3 is a view for explaining the finite element method; and  
         [0025]    [0025]FIG. 4 is a view for explaining the finite element method. 
     
    
     DETAILED DESCRIPTION OF THE PREFERRED EMBODIMENTS  
       [0026]    Preferred embodiments of the present invention will now be described in detail in accordance with the accompanying drawings.  
         [0027]    &lt;Apparatus Arrangement&gt; 
         [0028]    [0028]FIG. 1 is a block diagram showing an example of a processing apparatus which comprises a finite element method library according to an embodiment of the present invention. A programmer programs a finite element method program using this processing apparatus.  
         [0029]    Referring to FIG. 1, reference numeral  101  denotes a display device such as a CRT or the like; and  102 , various input devices such as a keyboard, mouse, or the like. Various programs are stored in a hard disk  103 . Such programs are loaded onto a memory  105  via a bus  107 , and are processed by a CPU  106 .  
         [0030]    The programs stored in the hard disk  103  are required to generate a finite element method program. Reference numeral  103 - 6  denotes an OS;  103 - 2 , an editor used to describe a program;  103 - 4 , various libraries;  103 - 5 , a header file;  103 - 1 , a compiler; and  103 - 3 , a linker.  
         [0031]    Note that these programs stored in the hard disk  103  may be stored in a floppy disk  108 , and may be loaded via a floppy disk drive  104 .  
         [0032]    &lt;Flow of Program Generation&gt; 
         [0033]    [0033]FIG. 2 shows the sequence in which the programmer programs a finite element method program using the above processing apparatus.  
         [0034]    Reference numeral  201  denotes a finite element method program (source) programmed using the editor  103 - 2 , and the finite element method program  201  has a header field  202  and a field  203  that describes the processing contents. The programmed finite element method program  201  is compiled using the compiler  103 - 2  to generate an object file  205 , to which the libraries  103 - 4  are linked by the linker  103 - 3 , thus generating an executable file  209 . In FIG. 2, reference numeral  206  denotes a finite element method library according to this embodiment of the plurality of libraries  103 - 4  to be linked by the linker  103 - 3 . More specifically, the finite element method library  206  comprises a header file  207  present as a text file, and a binary file  208  of the finite element method library. Of these files, the binary file  208  of the finite element method library is linked by the linker  103 - 3 , and the header file  207  is loaded upon being compiled by the compiler  103 - 2 .  
         [0035]    The finite element method library  206  is loaded from the hard disk  103  onto the memory  105  when it is used, in response to a signal from the input device  102  such as a keyboard, mouse, or the like input by the user, who develops a computer simulation program (finite element method program  201 ) which exploits the finite element method to analyze a physical phenomenon. As for the use method, the program is used by preparing its specification sheet of this embodiment. By launching the editor  103 - 2  using inheritance of classes in accordance with the specifications of respective classes, a program is described as an ASCII file, and the header file  207  (file name: FEM.h) upon compiling and the source  208  (file name: libFEM.a) of the finite element method library upon linking in generation of an execution program are loaded on the memory  105  together with the program when that program is used.  
         [0036]    With such flow of program generation, the programmer debugs bugs in the program on the basis of errors upon compiling or linking.  
         [0037]    &lt;Finite Element Method Library&gt; 
         [0038]    Details of the finite element method library  206  in FIG. 2 according to the embodiment of the present invention will be described below.  
         [0039]    In the finite element method library according to the embodiment of the present invention, as a type that represents a vector space, a type named FEMVector is defined in the header file  207  using a structure. Likewise, FEMCovector is prepared as vectors of a dual space. The present invention does not depend on individual type-declared names, but specific names are used in the above and following descriptions for the sake of simplicity.  
         [0040]    FEMGeometry is defined as a structure which has discrete geometric information that expresses a system. For example, structure FEMGeometry holds a sequence of FEMNode as structures which express lattice points.  
         [0041]    A function f(x) ( 204  in FIG. 2) described by the programmer in accordance with each individual problem is expressed by using basis e[N](x) that assume nonzero values in neighboring cells of each lattice point N and zero in other cells, and distributing real numbers f[N] to the basis. That is, f(x) is expressed by the linear combination given by:  
               f        (   x   )       =       ∑     all                 lattice                 points                 N              f        [   N   ]            e        [   N   ]            (   x   )                 (   4   )                               
 
         [0042]    A set {f[N]| N is all lattice points} of these real numbers distributed to the basis is considered as elements f=(f[N]) of the vector space V in this specification. FEMVector is defined as a structure that contains this set {f[N]| N is all lattice points} as a string. In this manner, the function f and vector space {f[N]| N is all lattice points} defined on the system, and FEMVector as the contents of the memory  105  on the processing apparatus can be identified. Since there is no confusion about symbols, it is expressed by f.  
         [0043]    On the other hand, a natural inner product  
               M        [     N   ,     N   ′       ]       =       ∫     entire                 domain                 of                 system                                           n            xe        [   N   ]              (   x   )          e        [     N   ′     ]            (   x   )                 (   5   )                               
 
         [0044]    is defined between neighboring basis, and it is natural to consider this inner product as a metric as linear mapping from the vector space V to the dual vector space V* (N′ indicates a lattice point different from N). Note that n is a quantity that represents a dimension which assumes 3 if the system is a three-dimensional system or 2 if the system is a two-dimensional system. Due to the linearity of the space, a matrix M is determined by only the natures of geometric information and basis e[N], and does not depend on individual vectors f. Hence, a dual vector f* (=f*[N]) of the vector f is defined by:  
                 f   *          [   N   ]       =         ∑     all                 lattice                 points                   N   ′                M        [     N   ,     N   ′       ]            f        [     N   ′     ]           =       ∫     entire                 domain                 of                 system                                           n            xe        [   N   ]              (   x   )          f        (   x   )                     (   6   )                               
 
         [0045]    In the finite element method library  206  of this embodiment, a new structure named FEMCovector is defined on the definition field  207  to form the quantity defined as described above on the memory  105  of the computer. That is, the matrix M is implemented as means for converting FEMVector into FEMCovector. That is, in the finite element method library  206  of this embodiment,  
                 DualMap        (   f   )       =     f   *       ,         f   *          [   N   ]       =       ∑     all                 lattice                 points                 N              M        [     N   ,     N   ′       ]            f        [     N   ′     ]                     (   7   )                               
 
         [0046]    is defined as a function. In general, since the calculation of an inverse matrix of M requires high calculation cost, no inverse operation of the above operation is prepared in the finite element method library  206  of this embodiment. However, such function may be prepared, or a function that implements the inverse operation by approximation may be prepared.  
         [0047]    As already described, an operation of a desired partial differential operator L to the function space F is limited to the Hilbert space H used by the finite element method library  206  of this embodiment, and is defined via:  
                 (   L   )          [     N   ,     N   ′       ]       =       ∫     entire                 domain                 of                 system                                           n            xe        [   N   ]              (   x   )          Le        [     N   ′     ]            (   x   )                 (   8   )                               
 
         [0048]    In consideration of this fact, the operation of the desired partial differential operator L to the function f of the system of this embodiment is implemented in the program as:  
                 (   Lf   )          [   N   ]       =       ∫     entire                 domain                 of                 system                                           n            xe        [   N   ]              (   x   )          Lf        (   x   )                   (   9   )                               
 
         [0049]    That is, an operation of L to f is defined as a conversion matrix from FEMVector to FEMCovector, and an operation of L to f is stored as FEMCovector on the memory  105 .  
         [0050]    In the finite element method library of this embodiment, pairing of such natural vector space and its dual space is included as a global function Pairing(FEMVector f, FEMCovector G). That is, the global function Pairing(FEMVector f, FEMCovector G) is defined by:  
               Pairing        (     f   ,   G     )       =       ∑     all                 lattice                 points                 N              f        [   N   ]            G        [   N   ]                   (   10   )                               
 
         [0051]    If elements of the vector space V are defined to have G =g* (i.e., g=M-1 G), a natural inner product  
               Pairing        (     f   ,     g   *       )       =       〈     f   ,   g     〉     =       ∫     entire                 domain                 of                 system                                           n            xf        (   x   )              g        (   x   )                     (   11   )                               
 
         [0052]    derived from the square integral (L2) norm of the function space (Iwanami&#39;s Dictionary of Mathematics, third edition, Iwanami, p. 61, 1994) can be implemented using f and g as the functions of the system. Hence, the finite element method library  206  of this embodiment provides natural pairing “Pairing” alone in place of the inner products of vector spaces or those of dual spaces.  
         [0053]    Using this Pairing, an inner product like &lt;f,Lg&gt; can be calculated. When (L) [N, N′] and M[N, N′] are uniquely determined if the geometry and basis function are determined, and the carrier of the basis function is locally defined, matrix elements can be calculated very easily. In practice, the domains of integrals that appear in M[N, N′] and (L) [N, N′] are on the entire domain. However, if U[N] is defined as a carrier of e[N] due to locality of the support (closure of a domain of definition having nonzero values) of the basis function, integrals of equations (5) and (8) are respectively changed to:  
                 (   L   )          [     N   ,     N   ′       ]       =       ∫       U        [   N   ]       ⋂     U        [     N   ′     ]                                               n            xe        [   N   ]              (   x   )          Le        [     N   ′     ]            (   x   )                 (   12   )                 M        [     N   ,     N   ′       ]       =       ∫       U        [   N   ]       ⋂     U        [     N   ′     ]                                               n            xe        [   N   ]              (   x   )          e        [     N   ′     ]            (   x   )                 (   13   )                               
 
         [0054]    That is, the integral can be calculated only in the intersection of supports, and can be simply expressed in general. Especially, implementation in the structured lattices or their hierarchical expression is very easy.  
         [0055]    As can be seen from equations (12) and (13), (L) [N, N′] and M[N, N′] are associated with the volume of U[N]∠U[N′]. Especially, it is apparent that M[N, N′] has positive correlation if a non-negative value function is used as e[N]. For this reason, in the finite element method in which the volumes of U[N]∠U[N′] are not constant, the ratio of each nonzero element of the metric M[N, N′] is normally very large, and approximation of it to a unit matrix considerably impairs the precision, convergence speed, and the like of the calculation results.  
         [0056]    [First Embodiment] 
         [0057]    In the first embodiment, an example of the header file  207  of the finite element method library  206  described in the above embodiment will be explained.  
         [0058]    The finite element method library of this embodiment is configured based on C++ complying with the standards of ISO and ANCI in 2000. As a matter of course, the configuration of the above embodiment can obtain the same result if another programming language is used.  
         [0059]    For this purpose, as a type that expresses a vector space, a type FEMVector is defined using an object-oriented class. Likewise, FEMCovector is prepared as a vector of a dual space.  
         [0060]    Parts in the source of the finite element method library of this embodiment, which are associated with the above embodiment, are listed below.  
         [0061]    As a class having lattice point information, FEMNode is defined by:  
                                                                                                                                               Class FEMNode{                double   x, y, z;           int   i:           int   njs;           int*   j;           FEMNode*   Next;            public:                FENNode( double x, double y, double z,                int I, int j1, int j2, int j3,           int j4, int j5=−1, int j6=−1 );                int geti( ){return i;}           . . .            };                FEMGeometry that totalizes them is prepared.            class FEMGeometry{                int   Nsize;           int   Csize;           FEMNode*   Np;           FEMCell*   Cp;            public:                . . .            };                  
 
         [0062]    Using these classes, a class that includes a real number sequence is defined as FEMVector.  
                                                                                                               class FEMVector{                FEMGeometry*   Geomp;           double*   Vp;            public:                FEMVector(FEMGeotmetry&amp; Geom, const int&amp;           Nodenumber) {                Geomp = &amp;Geom;           Vp = new double [Nodenumber];                }           double&amp; operator( ) (FEMNode N){                return Vp[N.geti( )];                }           . . .            };                  
 
         [0063]    As described in the embodiment of the present invention, FEMCovector is defined by:  
                                                                                                                     class FEMCovector{                FEMGeometry*   Geomp;           double*   Vp;            public:                FEMCovector(FEMGeometry&amp; Geom, const int&amp;            Nodenumber) {                Geomp = &amp;Geom;           Vp = new double [Nodenumber];                }           double&amp; operator( ) (FEMNode N) {                return Vp[N.geti( )];                }           . . .            };                  
 
         [0064]    The following three global functions are prepared:  
         [0065]    FEMCovector DualMap(FEMVector&amp; f);  
         [0066]    FEMCovector Laplacian(FEMVector&amp; f);  
         [0067]    double Pairing(FEMVector&amp; f, FEMCovector&amp; g);  
         [0068]    The contents of these functions are respectively given by equations (7), (9) and (10). Note that Laplacian uses a Laplacian operator as a linear partial differential operator. Furthermore, these functions are given as global functions, but may be member functions of appropriate classes or may be capsulated using namespace or the like.  
         [0069]    In the finite element method library generated in this manner, the header file  206  is linked upon compiling, the binary field (the source  208  of the finite element method library) is linked upon linking and execution like in a normal library, and the finite element method library is used as a finite element method program used to solve various physical problems.  
         [0070]    A physical quantity inherited from FEMVector, e.g., temperature in case of an equation of heat conduction, is defined, and the equation is solved according to physical requests. For example, we can write:  
         [0071]    class Temperature :public FEMVector{ 
         [0072]    double unit;  
         [0073]    . . .  
         [0074]    } 
         [0075]    If the user who uses a library for the above type definition forms Pairing of FEMVector and FEMvector or FEMCovector and FEMCovector by mistake, a compile error occurs, and the user must recognize the metric space in the finite element method, thus avoiding precision drop and deficiency of reliability due to careless programming. For this reason, the reliability of the library for the calculation result itself can also be improved.  
         [0076]    [Second Embodiment] 
         [0077]    In the first embodiment, FEMVector and FEMCovector are considered as types, but they may be distinguished using namespace or a truth quantity. That is,  
         [0078]    enum SpaceType{Vector, Covector};  
         [0079]    is defined to define following FEMSpace in place of FEMVector and FEMCovector of the first embodiment:  
                                                                                                                     class FEMSpace{                FEMGeometry*   Geomp;           double*   Vp;           SpaceType   Dual            public:                FEMspace(FEMGeometry&amp; Geom, const int&amp;            Nodenumber,       SpaceType&amp; D) {                Geomp = &amp;Geom;           VP = new double [Nodenumber];           Dual = D;                }           double&amp; operator( ) (FEMNode N) {                return Vp[N.geti( )];                }            };                  
 
         [0080]    In correspondence with this,  
         [0081]    FEMSpace DualMap(FEMSpace&amp; f);  
         [0082]    FEMSpace Laplacian(FEMSpace&amp; f);  
         [0083]    double Pairing(FEMSpace&amp; f, FEMSpace&amp; g);  
         [0084]    are given as global functions. Note that DualMap(f and Laplacian(f) execute condition processes for, e.g., generating an error and terminating a program, and so forth if SpaceType as a member of their argument f is not Vector, and set SpaceType of FEMSpace as a return value to be Covector if f is Vector. Also, Pairing(f, g) has a conditional statement that a program does not run upon execution if SpaceType of f and g are Vector and Vector; or Covector and Covector. The reason why DualMap (f) undergoes a measure for terminating a program if SpaceType of f is Covector is that it is generally difficult to calculate an inverse matrix of the metric M given by formula (1). Of course, DualMap may be defined as mapping from a dual space to a space by approximation or strictly.  
         [0085]    Furthermore, these functions are given as global functions, but may be member functions of appropriate classes or may be capsulated using namespace or the like.  
         [0086]    With this configuration, the user recognizes, as an error, his or her misapprehension about the metric space as the subject of the present invention upon execution of a program. Hence, the reliability of the finite element method library and a program using that library can be improved. That is, problems such as calculation errors, an increase in convergence time, and the like due to user&#39;s misapprehension can be avoided.  
         [0087]    [Third Embodiment] 
         [0088]    The basic principle of the finite element method library of the present invention can be formed independently of programming languages. For example, a program can be generated in Fortran90 and C. Also, a program that can obtain the same effects can be programmed using Pascal or the like. Furthermore, a similar program can be generated using an STL function of C++.  
         [0089]    [Fourth Embodiment] 
         [0090]    In the description of the present invention, a physical quantity is expressed by a real number value function. However, the present invention can be applied to a case wherein a physical quantity is expressed by a multi-dimensional vector value or complex number value. Especially, using a language that is compatible to generic programming in high-level languages such as C++ or the like, simple implementation can be made.  
         [0091]    As described above, according to the present embodiment, the reliability of a program that uses the finite element method as a library can be improved, and calculation errors, an increase in convergence time, and the like due to programming errors can be avoided.  
         [0092]    As many apparently widely different embodiments of the present invention can be made without departing from the spirit and scope thereof, it is to be understood that the invention is not limited to the specific embodiments thereof except as defined in the claims.