femdisk.cpp

FreeFEM is an implementation of the GFEM language dedicated to the finite element method. It provid

// Emacs will be in -*- Mode: c++ -*-
//
// SUMMARY: Language for a Finite Element Method
//
// ORG    : SOURCEFORGE
// E-MAIL : [e-mail]
//
// ORIG-DATE:     June-94
// LAST-MOD: 27-Oct-01 at 12:14:23 by Christophe Prud'homme
//
// DESCRIPTION:

/*
  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; either version 2 of the License, or
  (at your option) any later version.
  
  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., 675 Mass Ave, Cambridge, MA 02139, USA.

*/

// DESCRIP-END.
//
#if defined(__GNUG__)
#pragma implementation
#endif /* __GNUG__ */

#include <cstdio>
#include <cstdlib>
#include <cstring>

#include <iostream>
#include <fstream>


#include <femMisc.hpp>
#include <femDisk.hpp> 

#define NEXTLINE(ETIQ) {char c;while (((c = ETIQ.get()) != '\n') && !ETIQ.eof());}
#define END_OF_LINE(data,i,incr) { if (dummy == incr) {dummy = 0;fprintf(data,"\n");} else dummy ++;}

namespace fem
{
  extern int N;

  char *
  readprog (char *path)
  {
    long count;
    int l = 1;
    FILE *f;
    char *tmp = NULL;

    if ((f = fopen (path, "r")) == NULL)
      {
        fprintf (stderr, "Freefem::readprog error : Cannot read %s\n", path);
        exit(-1);
      }
    count = 0;
    while (!feof (f))
      {
        fgetc (f);
        count++;
      }
    rewind (f);
    tmp = new char[count + 255];
    memset(tmp, 0, (count + 255)*sizeof(char));
    tmp[0] = '{';			/* place program between braces */
    while (!feof (f))
      {
        fgets (tmp + l, 255, f);
        /*fputs (tmp + l, stdout);*/
        l = strlen (tmp);
      }
    tmp[l] = '}';
    l++;
    fclose (f);

    return tmp;
  }

  
  int 
  loadtriangulation (femMesh *t, char *path)
  {
    int i;
    char *result = 0;
    int dummy;

    std::ifstream fin( path );

    if ( fin.fail() )
      {
        return -1;
      }
  
    /*
     * strstr is found in string.h
     *  => portability problems
     *
     * amdba format 
     */
    if (( result = strstr (path,".amdba")) != 0 ) 
      {
        int __nPoints;
        int __nTriangles;

        fin >> __nPoints >> __nTriangles;
        NEXTLINE (fin);

        t->Delete();
        t->setDimensions( __nPoints, __nTriangles );
      
        for (i = 0; i < t->getNumberOfPoints(); i++)
          {
            fin >> dummy >> t->rp[i][0] >> t->rp[i][1] >> t->ng[i];
          }
      
          

        for (i = 0; i < t->getNumberOfCells(); i++)
          {
            fin >> dummy >> t->tr[i][0] >> t->tr[i][1] >> t->tr[i][2] >> t->ngt[i];
            t->tr[i][0]--;
            t->tr[i][1]--;
            t->tr[i][2]--;
          }
      }
    /* am_fmt format */
    else if ((result = strstr (path,".am_fmt")) != NULL)
      {
        int __nPoints;
        int __nTriangles;
      
        fin >> __nPoints >> __nTriangles;
        NEXTLINE (fin);

        t->Delete();
        t->setDimensions( __nPoints, __nTriangles );

        for (i = 0; i < t->getNumberOfCells(); i++)
          {
            fin >> t->tr[i][0] >> t->tr[i][1] >> t->tr[i][2];
            t->tr[i][0]--;
            t->tr[i][1]--;
            t->tr[i][2]--;
          }
        for (i = 0; i < t->getNumberOfPoints(); i++)
          {
            fin >> t->rp[i][0] >> t->rp[i][1];
          }
        for (i = 0; i < t->getNumberOfCells(); i++)
          {
            fin >> t->ngt[i];
          }
        for (i = 0; i < t->getNumberOfPoints(); i++)
          {
            fin >> t->ng[i];
          }
      }
    else /* gfem format */
      {
        int __nPoints;
        int __nTriangles;
      
        fin >> __nPoints >> __nTriangles;
        NEXTLINE (fin);

        t->Delete();
        t->setDimensions( __nPoints, __nTriangles );
      
        for (i = 0; i < t->getNumberOfPoints(); i++)
          {
            fin >> t->rp[i][0] >> t->rp[i][1] >> t->ng[i];
          }
        for (i = 0; i < t->getNumberOfCells(); i++)
          {
            fin  >> t->tr[i][0] >> t->tr[i][1] >> t->tr[i][2] >> t->ngt[i];
            t->tr[i][0]--;
            t->tr[i][1]--;
            t->tr[i][2]--;
          }
      }
    return 0;
  }

  int 
  savetriangulation (femMesh *t, char *path)
  {
    int i;
    char *result = NULL;

    std::ofstream fout( path );

    if ( fout.fail() )
      {
        return 1;
      }
    /*
     * strstr is found in string.h
     *  => portability problems
     */
    if ((result = strstr (path,".amdba")) != NULL) /* amdba format */
      {
        fout << t->getNumberOfPoints() << " " << t->getNumberOfCells() << std::endl;

        for (i = 0; i < t->getNumberOfPoints(); i++)
          {
            fout << i+1 << " " << t->rp[i][0] << " " << t->rp[i][1] << " " << t->ng[i] << std::endl;
          }
      
      
        for (i = 0; i < t->getNumberOfCells(); i++)
          {
            fout << i+1
                 << " " << t->tr[i][0]+1
                 << " " << t->tr[i][1]+1
                 << " " << t->tr[i][2]+1
                 << " " << t->ngt[i] << std::endl;
          }
      }
    else if ((result = strstr (path,".am_fmt")) != NULL)/* am_fmt format */
      {
#if 0
        fout << t->getNumberOfPoints() << " " << t->getNumberOfCells() < std::endl;
      
        dummy = 0;
        for (i = 0; i < t->getNumberOfCells(); i++)
          {
            fout << t->tr[i][0]+1
                 << " " << t->tr[i][1]+1
                 << " " << t->tr[i][2]+1
                 << std::endl;
            END_OF_LINE(data, dummy, 1);
          }
        if (dummy) fprintf(data,"\n");
        dummy = 0;
        for (i = 0; i < t->getNumberOfPoints(); i++)
          {
            fout << t->rp[i][0] << " " << t->rp[i][1] << std::endl;
            END_OF_LINE (data, dummy, 1);
          }
        if (dummy) fprintf(data,"\n");
        dummy = 0;
        for (i = 0; i < t->getNumberOfCells(); i++)
          {
            fout << t->ngt[i] << std::endl;
            END_OF_LINE (data, dummy, 9);
          }
        if (dummy) fprintf(data,"\n");
        dummy = 0;
        for (i = 0; i < t->getNumberOfPoints(); i++)
          {
            fout << t->ng[i] << std::endl;
            END_OF_LINE (data, dummy, 9);
          }
#endif
      }
    else
      {
        fout << t->getNumberOfPoints() << " " << t->getNumberOfCells() << std::endl;

        for (i = 0; i < t->getNumberOfPoints(); i++)
          {
            fout << t->rp[i][0] << " " << t->rp[i][1] << " " << t->ng[i] << std::endl;
          }
      
      
        for (i = 0; i < t->getNumberOfCells(); i++)
          {
            fout << t->tr[i][0]+1
                 << " " << t->tr[i][1]+1
                 << " " << t->tr[i][2]+1
                 << " " << t->ngt[i] << std::endl;
          }
      }
    return 0;
  }
  int
  savefct (creal * f, int ns, char *path)
  {
    char *result = NULL;
    int i; 

    std::ofstream fout( path );
    if ( fout.fail() )
      return 1;

    fout.precision( 8 );
    if ((result = strstr (path,".bb")) != NULL)/* bb format */
      {
        fout << "3 1 " << ns << " 2\n";
        for (i = 0;i < ns;i++)
          {
            fout << realpart( f[i] ) << "\n";
          }
      }
    else /* gfem format */
      {
        fout << ns << "\n";
        for (i = 0;i < ns;i++)
          {
            fout << realpart( f[i] ) << "\n";
          }
      }
    return 0;
  }

  int 
  saveparam (fcts * param, femMesh * t, char *path, int N)
  {
    int k, ns = t->getNumberOfPoints();
    
    std::ofstream file (path);
    file.precision (8);
    file << ns << "	" << N << std::endl;
    for (k = 0; k < ns; k++)
      {
        if (N == 1)
          {

            file << (param)->f1[k] << " ";
            file << "		";
            file << (param)->g1[k] << " ";
            file << "		";
            file << (param)->p1[k] << " ";
            file << "		";
            file << (param)->b1[k] << " ";
            file << "		";
            file << (param)->c1[k] << " ";
            file << "		";
            file << (param)->a11[k] << " ";
            file << "		";
            file << (param)->a21[k] << " ";
            file << "		";
            file << (param)->nuxx1[k] << " ";
            file << "		";
            file << (param)->nuxy1[k] << " ";
            file << "		";
            file << (param)->nuyx1[k] << " ";
            file << "		";
            file << (param)->nuyy1[k] << " ";
            file << "		";
          }
        else if (N == 2)
          {
            file << (param)->f2[k] << " ";
            file << "		";
            file << (param)->g2[k] << " ";
            file << "		";
            file << (param)->p2[k] << " ";
            file << "		";
            file << (param)->b2[k] << " ";
            file << "		";
            file << (param)->c2[k] << " ";
            file << "		";
            file << (param)->a12[k] << " ";
            file << "		";
            file << (param)->a22[k] << " ";
            file << "		";
            file << (param)->nuxx2[k] << " ";
            file << "		";
            file << (param)->nuxy2[k] << " ";
            file << "		";
            file << (param)->nuyx2[k] << " ";
            file << "		";
            file << (param)->nuyy2[k] << " ";
            file << "		";
          }
        file << std::endl;
      }
    file.close ();
    return 0;
  }


  int 
  saveconst (creal f, char *path)
  {
    std::ofstream file (path, std::ios::out | std::ios::app);
    //  file.seekoff(0,ios::end,0);
    file << f << std::endl;
    file.close ();
    return 0;
  }

  int 
  loadfct (creal * f, int ns, char *path)
  {
    char *result = NULL;
    int n,i,dummy1,dummy2,dummy3;

    std::ifstream fin(path);
    if ( fin.fail() )
      return 0;
    
    if ((result = strstr (path,".bb")) != NULL)/* bb format */
      {
        //fscanf (data,"%d %d %d %d",&dummy1,&dummy2,&n,&dummy3);
        fin >> dummy1 >> dummy2 >> n >> dummy3;
        NEXTLINE(fin);
        if (n != ns)
          return 0;
        for (i = 0;i < ns;i++)
          {
            //fscanf (data,"%g",&(f[i].real()));
            fin >> f[i].real();
            NEXTLINE(fin);
          }
      }
    else /* gfem format */
      {
        //fscanf (data,"%d",&n);
        fin >> n;
        NEXTLINE(fin);
        if (n != ns)
          return 0;
        for (i = 0;i < ns;i++) 
          {
            //fscanf (data,"%g",&(f[i].real()));
            fin >> f[i].real();
            NEXTLINE(fin);
          }
      }
    return -2;
  }

  int readpoints(char *path, float* cr, int imax)
  {
    int i = 0;
    std::ifstream fin( path );

    if ( fin.fail() )
      return 0;
    while ( !fin.eof() && i<imax)
      {
        //fscanf (data,"%8f\t%8f",&(cr[2*i]), &(cr[2*i+1]));
        fin >> cr[2*i] >> cr[2*i+1];
        NEXTLINE(fin); i++;
      }
    
    if(i==imax)
      return -1;
    else
      return --i;
    return 0;
  }

}