creep_dpl.cpp

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

   TKr, 3.6.2005
*/
void dplmat::give_strains(vector &eps,long ipp,long ido)
{
  long i;
  
  nc=Mm->ip[ipp].ncompstr;
  
  for (i=0;i<nc;i++){
    eps[i] = Mm->ip[ipp].eqother[ido+(0)+i];
  }
}



/**
   function stores total strains into ip

   @param eps - vector of total strains
   @param nc       - number of components
   @param ipp      - number of integration point

   TKr, 3.6.2005
*/
void dplmat::store_strains(vector eps,long ipp,long ido)
{
  long i;
  
  nc=Mm->ip[ipp].ncompstr;
  
  for (i=0;i<nc;i++){
    Mm->ip[ipp].eqother[ido+(0)+i] = eps[i];
  }
}


/**
   function returns increments of total strains from ip

   @param deps - vector of increments of total strains
   @param nc   - number of components
   @param ipp  - number of integration point

   TKr, 3.6.2005
*/
void dplmat::give_dstrains(vector &deps,long ipp,long ido)
{
  long i;
  
  nc=Mm->ip[ipp].ncompstr;
  
  for (i=0;i<nc;i++){
    deps[i] = Mm->ip[ipp].eqother[ido+(nc+nc)+i];
  }
}



/**
   function stores increments of total strains into ip

   @param deps - vector of increments of total strains
   @param nc   - number of components
   @param ipp  - number of integration point

   TKr, 3.6.2005
*/
void dplmat::store_dstrains(vector deps,long ipp,long ido)
{
  long i;
  
  nc=Mm->ip[ipp].ncompstr;
  
  for (i=0;i<nc;i++){
    Mm->ip[ipp].eqother[ido+(nc+nc)+i] = deps[i];
  }
}

 
/**
   function returns increments of irreversible creep strains into ip

   @param deps_cr - vector of increments of strains
   @param nc      - number of components
   @param ipp     - number of integration point

   TKr, 3.6.2005
*/
void dplmat::give_creepdstrains(vector &deps_cr,long ipp,long ido)
{
  long i;
  long n_ret_times;
  
  nc=Mm->ip[ipp].ncompstr;
  n_ret_times = give_nret_time ();
  
  for (i=0;i<nc;i++){
    deps_cr[i] = Mm->ip[ipp].eqother[ido+(nc+nc+nc+nc+2+2+n_ret_times+nc*n_ret_times)+i];
  }
}


/**
   function stores increments of irreversible creep strains into ip

   @param deps_cr - vector of increments of strains
   @param nc      - number of components
   @param ipp     - number of integration point

   TKr, 3.6.2005
*/
void dplmat::store_creepdstrains(vector deps_cr,long ipp,long ido)
{
  long i;
  long n_ret_times;
  
  nc=Mm->ip[ipp].ncompstr;
  n_ret_times = give_nret_time ();
  
  for (i=0;i<nc;i++){
    Mm->ip[ipp].eqother[ido+(nc+nc+nc+nc+2+2+n_ret_times+nc*n_ret_times)+i] = deps_cr[i];
  }
}


/**
   function returns increments of irreversible shrinkage strains into ip

   @param deps_sh - vector of increments of strains
   @param nc      - number of components
   @param ipp     - number of integration point

   TKr, 3.6.2005
*/
void dplmat::give_irrdstrains(vector &deps_sh,long ipp,long ido)
{
  long i;
  long n_ret_times;
  
  nc=Mm->ip[ipp].ncompstr;
  n_ret_times = give_nret_time ();
  
  for (i=0;i<nc;i++){
    deps_sh[i] = Mm->ip[ipp].eqother[ido+(nc+nc+nc+nc+2+2+n_ret_times+nc*n_ret_times+nc)+i];
  }
}


/**
   function stores increments of irreversible shrinkage strains into ip

   @param deps_sh - vector of increments of strains
   @param nc      - number of components
   @param ipp     - number of integration point

   TKr, 3.6.2005
*/
void dplmat::store_irrdstrains(vector deps_sh,long ipp,long ido)
{
  long i;
  long n_ret_times;
  
  nc=Mm->ip[ipp].ncompstr;
  n_ret_times = give_nret_time ();
  
  for (i=0;i<nc;i++){
    Mm->ip[ipp].eqother[ido+(nc+nc+nc+nc+2+2+n_ret_times+nc*n_ret_times+nc)+i] = deps_sh[i];
  }
}

/**
   function stores actual humidity

   @param ipp     - number of integration point

   TKr, 3.6.2005
*/
void dplmat::store_hum(long ipp,long ido)
{
  nc=Mm->ip[ipp].ncompstr;
  
  Mm->ip[ipp].eqother[ido+(nc+nc+nc+nc)] = 0.0;
}



/**
   function stores actual temperature

   @param ipp     - number of integration point

   TKr, 3.6.2005
*/
void dplmat::store_temp(long ipp,long ido)
{
  nc=Mm->ip[ipp].ncompstr;
  
  Mm->ip[ipp].eqother[ido+(nc+nc+nc+nc+1)] = 0.0;
}


/**
   function returns time of concrete (end of concrete casting) in days

   @param ipp - number of integration point

   TKr, 3.6.2005
*/
double dplmat::give_tb_time(long ipp)
{
  double tb;
  
  tb = tb_time/86400.0;//days

  return(tb);
}



/**
   function returns time when drying begins in days

   @param ipp - number of integration point

   TKr, 3.6.2005
*/
double dplmat::give_th_time(long ipp)
{
  double th;
  
  th = th_time/86400.0;//days

  return(th);
}



/**
   function returns number of times in which complinace function is approximated

   @param ipp - number of integration point

   TKr, 3.6.2005
*/
long dplmat::give_napproxtime(long ipp)
{
  return(napproxtime);
}


void dplmat::nlstresses (long ipp,long im,long ido)
{
  long i,nm, ncompstr = Mm->ip[ipp].ncompstr, ncompo;
  vector sigback(ncompstr);

  if (Mp->phase==1)
  {
    // computation of increments of stresses due to the irreversible strain increment
    Mm->creep_phase1(ipp,im,ido);

    // computation of increments of stresses due to the temperature strain increment
    if ((im == 0) && (Mm->ip[ipp].hmt & 1))
    {
      nm=Mm->ip[ipp].nm-1;
      ncompo = Mm->givencompeqother(ipp, 0);
      ncompo -= Mm->givencompeqother(ipp, nm);
      for (i=0; i < ncompstr; i++)
        sigback[i] = Mm->ip[ipp].stress[i];
      Mm->computenlstresses(ipp, nm, ncompo);
      for (i=0; i < ncompstr; i++)
        Mm->ip[ipp].stress[i] += sigback[i];
    }
  }
  
  if (Mp->phase==2)
  {
    Mm->creep_phase2(ipp,im,ido);
    if ((im == 0) && (Mm->ip[ipp].hmt & 1))
    {
      nm=Mm->ip[ipp].nm-1;
      ncompo = Mm->givencompeqother(ipp, 0);
      ncompo -= Mm->givencompeqother(ipp, nm);
      // actualization of previous temperature
      Mm->computenlstresses(ipp, nm, ncompo);
    }
  }
}


/**
   function computes creep compliacne function J(t,t') [1/Pa] in a material point Double power law Bazant's model

   @param tt        - age at loading [days]
   @param t         - time representing age of concrete [days]
   @param ipp       - number of integration point
  
   TKr, 3.6.2005
*/
double dplmat::double_power_law(double tt,double t,long ipp)
{
  double ee,a,n,m,fi;
  double x,ag,ac,jt;
  
  if (type_e == 1){//measured
    ee = 1.0/e0*1.0e6;//ee must be in 1.0e-6/psi
  }
  else {
    ee=(0.09+1/(1.7*(0.00005*ro*ro*fc/1000.0)*(0.00005*ro*ro*fc/1000.0)));//fc must be in ksi
  }
  
  m=0.28+1/(fc*fc/1000.0/1000.0);
  a=1./(40.0*wc); 
  ag=(sc+gc)/gc;
  ac=sc+gc;
  x=(2.1*ac/pow(sc,1.4)+0.1*pow((fc/1000.0),1.5)*pow(wc,1./3.)*pow(ag,2.2))*a-4.0;//fc must be in ksi
  if (x <= 4.0)
    n=0.12;
  else
    n=.12+(0.07*pow(x,6.)/(5130.0+pow(x,6.)));
  
  fi=0.5*pow(10.0,(3.0*n))/(pow(28.0,-m)+a);
  
  jt=ee+fi*ee*(pow(tt,-m)+a)*pow((t-tt),n);

  jt = jt/6.89476*1.0e-3*1.0e-6;//units 1e-6*psi^-1 -> Pa^-1

  return(jt);
}


/**
   function computes time increment of irreversible free strains in a material point

   @param ipp       - number of integration point

   TKr, 3.6.2005
*/
void dplmat::give_deps_free (long ipp)
{
}



/**
   function computes time increment of irreversible stress induced strains in a material point
   
   @param ipp       - number of integration point

   TKr, 3.6.2005
*/
void dplmat::give_deps_stressinduced (long ipp)
{
}



/**
  Function stores irreversible strains.

  @param ipp   - integration point number in the mechmat ip array.
  @param ido   - index of the first internal variable for given material in the ipp other array
  @param deps  - %vector of increments of irreversible strains which to be stored into total irr. strains
 
  Returns vector of irreversible strains via parameter epscr
*/
void dplmat::storeirrstrains (vector deps,long ipp, long ido)
{
  long i;
  long n_ret_times;
  long nc;

  nc=Mm->ip[ipp].ncompstr;
  
  n_ret_times = give_nret_time ();
  fprintf (stdout,"\nstore deps = %e\n",deps[0]);
  for (i=0;i<deps.n;i++){
    Mm->ip[ipp].eqother[ido+(nc+nc+nc+nc+2+2+n_ret_times+nc*n_ret_times+nc+nc)+i] = Mm->ip[ipp].eqother[ido+(nc+nc+nc+nc+2+2+n_ret_times+nc*n_ret_times+nc+nc)+i] + deps[i];
  }
}



/**
  Function returns irreversible strains.

  @param ipp   - integration point number in the mechmat ip array.
  @param ido   - index of the first internal variable for given material in the ipp other array
  @param epscr  - %vector of irreversible strains
 
  Returns vector of irreversible strains via parameter epscr
*/
void dplmat::giveirrstrains (long ipp, long ido, vector &epscr)
{
  long i;
  long n_ret_times;
  long nc;
  
  nc=Mm->ip[ipp].ncompstr;
  
  n_ret_times = give_nret_time ();
  
  for (i=0;i<epscr.n;i++){
    epscr[i] = Mm->ip[ipp].eqother[ido+(nc+nc+nc+nc+2+2+n_ret_times+nc*n_ret_times+nc+nc)+i];
  }
  fprintf (stdout,"\n give deps = %e\n",epscr[0]);
}

/**
   function returns other value from ip
   @param compother - number of other components
   @param ipp - first integration point on element
*/

double dplmat::get_othervalue(long compother,long ipp)
{
  double other;

  //dodelat??!!
  switch (compother){
  case 0:{//eps_x
    other = 0.0;
      break;
  }
  default:{
    fprintf (stderr,"\n\n unknown type of component is required in function (%s, line %d).\n",__FILE__,__LINE__);
  }
  }
  return (other);

}

/**
     function prints name of other value
     @param out - output file
     @param compother - number of other components
*/
void dplmat::print_othervalue_name(FILE *out,long compother)
{
  //dodelat??!!
  switch (compother){
  case 0:{//eps_x
    fprintf (out,"eps_x ()");
    break;
  }
  default:{
    fprintf (stderr,"\n\n unknown type of component is required in function (%s, line %d).\n",__FILE__,__LINE__);
  }
  }
}