Line data    Source code

       1              : #ifndef MODEL_INTERFACE_CRACK_H
       2              : #define MODEL_INTERFACE_CRACK_H
       3              : 
       4              : #include <AMReX.H>
       5              : #include <AMReX_AmrCore.H>
       6              : 
       7              : #include <iostream>
       8              : #include <fstream>
       9              : 
      10              : namespace Model
      11              : {
      12              : namespace Interface
      13              : {
      14              : namespace Crack
      15              : {
      16              : class Crack
      17              : {
      18              :     public:
      19              :     enum GType {GSQUARE, GMULTIWELL, GPHI4C3, GCUBICM, GSQUAREP, GSQUAREPM};
      20              :     enum WType {WSQUARE, WMULTIWELL, WMULTIWELL2, WPHI4C3};
      21              :     
      22            0 :     Crack() {};
      23              :     
      24              :     AMREX_FORCE_INLINE
      25            0 :     virtual Set::Scalar w_phi(Set::Scalar c, Set::Scalar /*p=0.*/)
      26              :     {
      27            0 :         switch(w_type)
      28              :         {
      29            0 :             case WSQUARE:        return (1.-c)*(1.-c);
      30            0 :             case WMULTIWELL:        return (1.-c)*(1.-c)*c*c;
      31            0 :             case WMULTIWELL2:    return (1.+c)*(1.+c)*(1.-c)*(1.-c);
      32            0 :             case WPHI4C3:        return 1. - 4.*c*c*c + 3.*c*c*c*c;
      33              :             
      34            0 :             default:            return (1.-c)*(1.-c);
      35              :         }
      36              :     }
      37              :     AMREX_FORCE_INLINE
      38            0 :     virtual Set::Scalar g_phi(Set::Scalar c, Set::Scalar p=0.)
      39              :     {
      40            0 :         switch(g_type)
      41              :         {
      42            0 :             case GSQUARE:        return c*c;
      43            0 :             case GMULTIWELL:     return (2.-c)*(2.-c)*c*c;
      44            0 :             case GPHI4C3:         return 4.*c*c*c - 3.*c*c*c*c;
      45            0 :             case GSQUAREP:        return std::pow(c,2.*p);
      46            0 :             case GSQUAREPM:        return std::pow(c,2.*(std::pow(p,m_d_exponent)));
      47            0 :             case GCUBICM:        return m_d_exponent*(c*c*c - c*c) + 3.*c*c - 2.*c*c*c;
      48              : 
      49            0 :             default:             return c*c;
      50              :         }
      51              :     }
      52              :     AMREX_FORCE_INLINE
      53            0 :     virtual Set::Scalar Dw_phi(Set::Scalar c, Set::Scalar /*p=0.*/)
      54              :     {
      55            0 :         switch(w_type)
      56              :         {
      57            0 :             case WSQUARE:        return -2.*(1.-c);
      58            0 :             case WMULTIWELL:        return 4.*c*c*c - 6.*c*c + 2.*c;
      59            0 :             case WMULTIWELL2:    return 4.*c*c*c - 4.*c;
      60            0 :             case WPHI4C3:        return 12.*(c-1.)*c*c;
      61              :             
      62            0 :             default:            return -2.*(1.-c);
      63              :         }
      64              :     }
      65              :     AMREX_FORCE_INLINE
      66            0 :     virtual Set::Scalar Dg_phi(Set::Scalar c, Set::Scalar p=0.)
      67              :     {
      68            0 :         switch(g_type)
      69              :         {
      70            0 :             case GSQUARE:        return 2.*c;
      71            0 :             case GMULTIWELL:     return 4.*c*c*c - 12.*c*c + 8.*c;
      72            0 :             case GPHI4C3:         return 12.*(1.-c)*c*c;
      73            0 :             case GSQUAREP:        return 2.*p*std::pow(c,2*p -1.);
      74            0 :             case GSQUAREPM:        return 2.*std::pow(p,m_d_exponent)*(std::pow(c, 2*std::pow(p,m_d_exponent)-1));
      75            0 :             case GCUBICM:        return m_d_exponent*(3.*c*c - 2.*c) + 6.*c - 6.*c*c;
      76              : 
      77            0 :             default:             return 2*c;
      78              :         }
      79              :     }
      80              : 
      81              :     virtual Set::Scalar Gc(Set::Scalar theta) = 0;
      82              :     virtual Set::Scalar DGc(Set::Scalar theta) = 0;
      83              :     virtual Set::Scalar DDGc(Set::Scalar theta) = 0;
      84              :     virtual Set::Scalar Zeta(Set::Scalar theta) = 0;
      85              :     virtual Set::Scalar Mobility(Set::Scalar theta) = 0;
      86              :     virtual Set::Scalar DrivingForceThreshold(Set::Scalar theta) = 0;
      87              : 
      88              :     void ExportToFile(std::string filename, amrex::Real dTheta)
      89              :     {
      90              :         std::ofstream outFile;
      91              :         outFile.open(filename);
      92              :     
      93              :         for(amrex::Real theta=0; theta<2*pi; theta=theta+dTheta)
      94              :             {
      95              :                 outFile << theta << " " << Gc(theta) << std::endl;   
      96              :             }
      97              :         outFile.close();
      98              :     
      99              :     }
     100            0 :     void SetGType(const GType a_type)
     101              :     {
     102            0 :         g_type = a_type;
     103            0 :     }
     104              : 
     105            0 :     void SetWType(const WType a_type)
     106              :     {
     107            0 :         w_type = a_type;
     108            0 :     }
     109              :     
     110            0 :     void SetDuctileExponent(const Set::Scalar m)
     111              :     {
     112            0 :         m_d_exponent = m;
     113            0 :     }
     114              : 
     115              : protected:
     116              :     static constexpr amrex::Real pi = 3.14159265359;
     117              :     GType g_type = GType::GSQUARE;
     118              :     WType w_type = WType::WSQUARE;
     119              :     Set::Scalar m_d_exponent = 1.;
     120              : };
     121              : }
     122              : }
     123              : }
     124              : 
     125              : #endif