main.cpp

来自「OFELI is an object oriented library of C」· C++ 代码 · 共 96 行

/*==============================================================================

                                    O  F  E  L  I

                            Object  Finite  Element  Library

  ==============================================================================

   Copyright (C) 1998 - 2004 Rachid Touzani

   This program is free software; you can redistribute it and/or modify it under
   the terms of the GNU General Public License as published by the Free 
   Software Foundation; Version 2 of the License.

   This program is distributed in the hope that it will be useful, but WITHOUT
   ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
   FOR A PARTICULAR PURPOSE. See the GNU General Public License for more 
   details.

   You should have received a copy of the GNU General Public License 
   along with this program; if not, write to the :

   Free Software Foundation
   Inc., 59 Temple Place - Suite 330
   Boston, MA  02111-1307, USA

  ==============================================================================

              An example of a Finite Element Code using OFELI

           Solution of a 2-D Diffusion problem using P1 Finite elements
           Linear System is solved using the Conjugate Gradient Method

  ==============================================================================*/

#include "OFELI.h"
#include "Therm.h"

using namespace OFELI;

int main(int argc, char *argv[])
{
   Element *el;
   banner();

// Read and output mesh data
   if (argc <= 1) {
     cout << " Usage : ex3 <mesh_file>" << endl;
     exit(1);
   }
   Mesh ms(argv[1]);
   cout << ms;
   ms.EliminateImposedDOF();
   ms.NumberEquations();

// Declare problem data (matrix, rhs, boundary conditions, body forces)
   SpMatrix<double> a(ms);
   Vect<double> b(ms.NbEq()), x(ms.NbEq());
   BCVect<double> bc(ms.BC());

// Loop over elements
// ------------------

   for (ms.TopElement(); (el=ms.GetElement());) {

//    Declare an instance of class DC2DT3
      DC2DT3 eq(el);

//    Diffusion contribution to matrix
      eq.Diffusion();

//    Boundary condition contribution to RHS
      eq.UpdateBC(bc);

//    Assemble element matrix and RHS
      a.Assembly(el,eq.A());
      b.Assembly(el,eq.b());
   }

// Solve the linear system of equations by an iterative method
   DiagPrec<double,SpMatrix<double> > p(a);
   double toler = 1.e-8;
   int nb_it = CG(a,p,b,x,1000,toler,0);
   cout << "Nb. of iterations : " << nb_it << endl;

// Output solution
   Vect<double> u(ms.NbDOF());
   u.InsertBC(ms,x,bc);
   cout << "\nSolution :\n" << u;
   
#ifdef WITH_PAUSE
   system("PAUSE");
#endif
   return 0;
}