damplast.cpp

Finite element program for mechanical problem. It can solve various problem in solid problem

#include "damplast.h"
#include "matrix.h"
#include "vector.h"
#include "elastisomat.h"
#include "global.h"
#include "intpoints.h"
#include "vecttens.h"
#include <math.h>

/**
  This constructor inializes attributes to zero values.
*/
damplast::damplast (void)
{

}

/**
  This destructor is only for the formal purposes.
*/
damplast::~damplast (void)
{

}

/**
  This function computes material stiffnes %matrix.

  @param d - allocated matrix structure for material stiffness %matrix
  @param ipp - integration point number

*/
void damplast::matstiff (matrix &d,long ipp)
{
  double dp;
  long ncompo = Mm->givencompother(ipp, 1);

  switch (Mp->stmat)
  {
    case initial_stiff :
      Mm->elmatstiff (d,ipp);
      break;
    case tangent_stiff :
      Mm->elmatstiff (d,ipp);
      dp=Mm->ip[ipp].eqother[ncompo+1];
      if (dp > 0.999999)
        dp = 0.999999;
      cmulm (1.0-dp,d);
      break;
    default :
      fprintf(stderr, "\n\nError - unknown type of stifness matrix");
      fprintf(stderr, "\n in function damplast::matstiff (%s, line %d)\n", __FILE__, __LINE__);
  }
}


/**
  This function computes correct stresses in the integration point and stores
  them into ip stress array.

  @param ipp - integration point pointer
  7.10.2001
*/
void damplast::nlstresses (long ipp)
{
  long i;
  long ncompo = Mm->givencompother(ipp, 1);
  long ncomp = Mm->ip[ipp].ncompstr;
  vector epsback(ncomp);
  vector kappa(1);
  vector sigma(ncomp);
  double dp;

  // compute plasticity, it has index 1
  Mm->computenlstresses(ipp, 1, 0);  
  
  //  backup of the total strains
  for (i=0;i<ncomp;i++){
    epsback[i] = Mm->ip[ipp].strain[i];
  }
  //  initial values of elastic strain for damage = total strain - plastic strain
  for (i=0;i<ncomp;i++){
    Mm->ip[ipp].strain[i] -= Mm->ip[ipp].other[i];
  }
  
  kappa[0] = Mm->ip[ipp].eqother[ncompo+0];

  // damage stress solver, damage material has index 2
  dp = Mm->scal_dam_sol (ipp, 2, ncompo, epsback, kappa, sigma);
  
  //  recovery of the total strains
  for (i=0;i<ncomp;i++){
    Mm->ip[ipp].strain[i] = epsback[i];
  }
  
  
  //  new data storage
  for (i=0;i<ncomp;i++)
    Mm->ip[ipp].stress[i]=sigma[i];
  
  Mm->ip[ipp].other[ncompo+0]=kappa[0];
  Mm->ip[ipp].other[ncompo+1]=dp;
  
}

/**
  This function updates values in the other array reached in the previous equlibrium state to
  values reached in the new actual equilibrium state.

  @param ipp - integration point number in the mechmat ip array.

*/
void damplast::updateval (long ipp)
{
  //  update values of plasticity model
  Mm->updateipvalmat (ipp,1,0);
  //  update values of damage model
  Mm->updateipvalmat (ipp,2,Mm->givencompother(ipp,1));
}