femmisc.hpp

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

    friend std::ostream& operator<<(std::ostream& os,  cvect& a);
    friend std::istream& operator>>(std::istream& os,  cvect& a);
  };

  //_______________________________
  class cmat { 
    //---------
  public: ccreal  val[4];       //  Internal cvector to n values
    cmat()  {val[0]=0.F;val[1]=0.F;val[2]=0.F;val[3]=0.F; }         // constructor
    cmat(const cmat& r) {val[0]=r.val[0];val[1]=r.val[1];val[2]=r.val[2];val[3]=r.val[3];}  // copy constructor
    ///cmat(cmat& r) { val[0]=r.val[0];  val[1]=r.val[1];  val[2]=r.val[2];  val[3]=r.val[3]; }  // copy constructor
    cmat(ccreal r) {  val[0] = r;val[1] = r;val[2] = r;val[3] = r;} // copy constructor from a creal produces a diag mat
    ccreal& operator() (int i, int j)   { /*myassert((i< N)&&(i>=0)&&(j<N)&&(j>=0));*/ return val[i+i+j];}
    ccreal& operator[] (int i)    { /*myassert((i< N2)&&(i>=0));*/  return val[i];}
    const ccreal& operator[] (int i) const    { /*myassert((i< N2)&&(i>=0));*/  return val[i];}
    float modul2()    {return val[0]+val[1]+val[2]+val[3];} // public fnct
    cmat operator*=(const cmat& a)  
    { cmat b;  
    b.val[0] =val[0]*a.val[0]+val[1]*a.val[2];
    b.val[1] =val[0]*a.val[1]+val[1]*a.val[3];
    b.val[2] =val[2]*a.val[0]+val[3]*a.val[2];
    b.val[3] =val[2]*a.val[1]+val[3]*a.val[3];
    return b; 
    }   
    cmat& operator*=(const ccreal a) { val[0] *= a; val[1] *= a; val[2] *= a; val[3] *= a; return *this; }   
    friend cmat operator *(const cmat& a, const cmat& b) { cmat c(a); return c *= b; }  
    friend cmat operator *(const cmat& a, const ccreal b) { cmat c(a); return c *= b; }  
    friend cmat operator *(const ccreal b, const cmat& a) { cmat c(a); return c *= b; }   
    cmat& operator/=(const ccreal a) { val[0] /= a; val[1] /= a; val[2] /= a; val[3] /= a; return *this; }   
    friend cmat operator /(const cmat& a, const ccreal b) { cmat c(a); return c /= b; }  
    cmat& operator+=(const cmat& a){val[0]+=a.val[0];val[1]+=a.val[1];val[2]+=a.val[2];val[3]+=a.val[3]; return *this; }   
    cmat& operator+=(const ccreal a) { val[0] += a; val[1] += a; val[2] += a; val[3] += a; return *this; }   
    friend cmat operator +(const cmat& a, const cmat& b) { cmat c(a); return c += b; }  
    friend cmat operator +(const cmat& a, const ccreal b) { cmat c(a); return c += b; }  
    friend cmat operator +(const ccreal b, const cmat& a) {  cmat c(a); return c += b; }   

    cmat& operator-=(const cmat& a){val[0]-=a.val[0];val[1]-=a.val[1];val[2]-=a.val[2];val[3]-=a.val[3]; return *this; }   
    cmat& operator-=(const ccreal a)  { val[0] -= a; val[2] -= a; val[1] -= a; val[3] -= a; return *this; }   
    friend cmat operator -(const cmat& a, const cmat& b) { cmat c(a); return c -= b; }  
    friend cmat operator -(const cmat& a, const ccreal b) { cmat c(a); return c -= b; }  
    friend cmat operator -(const ccreal b, const cmat& a) {  cmat c(b); return c -= a; }   
    cmat& operator=(const cmat& a){val[0]=a.val[0];val[1]=a.val[1];val[2]=a.val[2]; val[3] = a.val[3]; return *this; }     
    cmat& operator=(const ccreal a)   { val[0] = a; val[1] = a; val[2] = a; val[3] = a; return *this; }     
    cmat& operator-()   { val[0] = -val[0];val[1] = -val[1];val[2] = -val[2];val[3] = -val[3]; return *this; }     
    friend float norm2( cmat& a) { return a.modul2();}
    /**/
    friend cvect operator *(const cmat& a, const cvect& b) 
    {   cvect c;  c.val[0] = a.val[0]*b.val[0]+a.val[1]*b.val[1];
    c.val[1] = a.val[2]*b.val[0]+a.val[3]*b.val[1]; 
    return c ; }  
    friend cvect operator *(const cvect& b, const cmat& a)
    { cvect c;    c.val[0] = a.val[0]*b.val[0]+a.val[2]*b.val[1];
    c.val[1] = a.val[1]*b.val[0]+a.val[3]*b.val[1]; 
    return c ; }  //transposed
    friend cvect operator /(const  cvect& b,const cmat& a) 
    {       cvect c; ccreal s = deter2(a.val[0],a.val[1], a.val[2],a.val[3]);
    if(norm2(s)<epsilon) { s = epsilon;}
    c.val[0] =  deter2(b.val[0],a.val[1], b.val[1],a.val[3]) / s;
    c.val[1] =  deter2(a.val[0],b.val[0], a.val[2],b.val[1]) / s;
    return c;
    }  
    friend cmat operator /(const cmat& b,const cmat& a) 
    {
      cmat c; ccreal s = deter2(a.val[0],a.val[1], a.val[2],a.val[3]);
      if(norm2(s)<epsilon){ s = epsilon;}
      c.val[0] =  deter2(b.val[0],a.val[1], b.val[2],a.val[3]) / s;
      c.val[2] =  deter2(a.val[0],b.val[0], a.val[2],b.val[2]) / s;
      c.val[1] =  deter2(b.val[1],a.val[1], b.val[3],a.val[3]) / s;
      c.val[3] =  deter2(a.val[0],b.val[1], a.val[2],b.val[3]) / s;
      return c;
    } 
    friend std::ostream& operator<<(std::ostream& os,  cmat& a);
    friend std::istream& operator>>(std::istream& os,  cmat& a);
    friend cmat id(const cvect& a); 
  };

  /*
    class Afloat
    { public:     
    int szz; 
    float* cc; 
  
    Afloat():szz(0),cc(NULL){}
    Afloat(int sz=0)  {cc = 0;   
    if(sz>0){ cc = new float[sz]; 
    if(!cc) erreur("Out of Memory");
    for(int i=0;i<sz;i++)cc[i]=0; } 
    szz = sz;
    }
    Afloat(Afloat& a)       {   cc = 0;  szz = 0; 
    if(a.szz>0) { szz = a.szz; cc = new float[szz]; 
    if(!cc) erreur("Out of Memory");
    else for(int i=0;i<szz;i++)cc[i]=a.cc[i];} 
    }
    ~Afloat()           { delete  [] cc;cc=0;szz = 0;} 
    float& operator[] (int i)   { myassert((i< szz)&&(i>=0)); return cc[i];}
    void init(int newSize)    { myassert(!(szz || cc)); szz =newSize; 
    cc =new float[szz]; if(!cc) erreur("Out of Memory");
    for(int i=0;i<szz;i++)cc[i]=0;;   
    } 
    };
  */
  class Acreal
  { 
  public:   
    long szz; 
    creal* cc; 
  
    Acreal(long=0);
    Acreal(const Acreal&);
    ~Acreal()  { delete  [] cc;cc=0;szz = 0;} 
    void destroy()  {delete  [] cc;cc=0;szz = 0;} 
    creal& operator[] (long i)  { /*myassert((i< szz)&&(i>=0));*/ return cc[i];}
    creal* operator&()      { return cc;} 
    void init(long newSize);/*    { myassert(!(szz || cc)); szz =newSize; 
				  cc =new creal[szz];   if(!cc) erreur("Out of Memory");
				  for(int i=0;i<szz;i++)cc[i]=0; 
				  } */
  };
  class Aint
  { 
  public:
    long szz; 
    int* cc; 
  
    Aint(long=0);
    Aint(const Aint&);
    ~Aint()           { delete  [] cc;cc=0;szz = 0;}
    void destroy()  {delete  [] cc;cc=0;szz = 0;} 
    int& operator[] (long i)  { myassert((i< szz)&&(i>=0)); return cc[i];}
    int* operator&()      { return cc;} 
    void init (long);
  };


  class Acvect
  {
  public:
    long szz; 
    cvect* cc; 
  
    Acvect(long=0);
    Acvect(const Acvect&);
    ~Acvect()           { delete  [] cc;cc=0;szz = 0;} 
    void destroy()  { delete  [] cc;cc=0;szz = 0;} 
    cvect& operator[] (long i)  { /*myassert((i< szz)&&(i>=0));*/ return cc[i];}
    cvect* operator&()      { return cc;} 
    void init (long);
  };

  class Acmat
  { 
  public: 
    long szz; 
    cmat* cc; 
  
    Acmat(long=0);
    Acmat(const Acmat&);
    ~Acmat()          
    {
      delete  [] cc;
      cc=0;
      szz = 0;
    } 
    void destroy()  {delete  [] cc;cc=0;szz = 0;} 
    cmat& operator[] (long i)   { /*myassert((i< szz)&&(i>=0));*/ return cc[i];}
    cmat* operator&()     { return cc;} 
    void init (long);
  };

  class AAcmat
  {
  public:
    long szz; 
    Acmat* cc; 

    AAcmat(long=0);
    AAcmat(const AAcmat&);
    ~AAcmat()           { delete  [] cc;cc=0;szz = 0;} 
    void destroy()  {delete  [] cc;cc=0;szz = 0;} 
    Acmat& operator[] (long i)  { /*myassert((i< szz)&&(i>=0));*/ return cc[i];}
    Acmat* operator&()      { return cc;} 
    void init(long );
  };

  class AAcreal
  { 
  public:
    long szz; 
    Acreal* cc; 
  
    AAcreal(long=0);
    AAcreal(const AAcreal& );
    ~AAcreal()          { delete  [] cc;cc=0;szz = 0;} 
    void destroy()  {delete  [] cc;cc=0;szz = 0;} 
    Acreal& operator[] (long i)   { /*myassert((i< szz)&&(i>=0));*/ return cc[i];}
    Acreal* operator&()     { return cc;} 
    void init (long);
  };
}

#endif /* __OPCLASS_H */