00001 #include <math.h>
00002 #include "cebfipcontactmat.h"
00003 #include "matrix.h"
00004 #include "vector.h"
00005 #include "stochdriver.h"
00006 #include "global.h"
00007 #include "intpoints.h"
00008
00009 cebfipcontactmat::cebfipcontactmat (void)
00010 {
00011 fcc = 0.0; fct = 0.0; conf = 0; bond = 0; normal = 0.0;
00012
00013 s = 0.0; ss = 0.0; sss = 0.0; taumax = 0.0; tauf = 0.0; alfa = 0.0;
00014 }
00015
00016 cebfipcontactmat::~cebfipcontactmat (void)
00017 {
00018
00019 }
00020
00021 void cebfipcontactmat::read (XFILE *in)
00022 {
00023 xfscanf (in,"%lf %lf %ld %ld %lf",&fcc,&fct,&conf,&bond,&normal);
00024
00025 if (conf==1)
00026 {
00027 s=ss=0.6e-3; alfa=0.4;
00028 if (bond==1) {sss=1.0e-3;taumax=2*sqrt(fcc/1e6)*1e6;tauf=0.3*sqrt(fcc/1e6)*1e6;}
00029 if (bond==2) {sss=2.5e-3;taumax=sqrt(fcc/1e6)*1e6;tauf=0.15*sqrt(fcc/1e6)*1e6;}
00030 }
00031 if (conf==2)
00032 {
00033 s=1.0e-3; ss=3.0e-3; alfa=0.4;
00034 if (bond==1) {taumax=2.5*sqrt(fcc/1e6)*1e6;tauf=1.0*sqrt(fcc/1e6)*1e6;}
00035 if (bond==2) {taumax=1.25*sqrt(fcc/1e6)*1e6;tauf=0.5*sqrt(fcc/1e6)*1e6;}
00036 }
00037 }
00038
00039 void cebfipcontactmat::matstiff (matrix &d,long ipp)
00040 {
00041 long i, n = Mm->ip[ipp].ncompstr;
00042 vector eps(n);
00043
00044 for (i=0;i<n;i++)
00045 eps[i]=Mm->ip[ipp].strain[i];
00046
00047 eps[0]=fabs(eps[0]);
00048 if (eps[0]<0.05*s)
00049 eps[0] = 0.05*s;
00050 if (eps[0]<0.8*s)
00051 d[0][0]=taumax/pow(s,alfa)*alfa*pow(eps[0],(alfa-1));
00052 else
00053 d[0][0]=taumax/pow(s,alfa)*alfa*pow(0.8*s,(alfa-1));
00054
00055 d[0][0] = d[0][0]*0.0314;
00056 d[1][1] = normal;
00057 }
00058
00059 void cebfipcontactmat::nlstresses (long ipp, long im, long ido)
00060 {
00061 long i, n = Mm->ip[ipp].ncompstr;
00062 vector eps(n),sig(n),dpp(n),dpmax(n);
00063
00064
00065 for (i=0;i<n;i++)
00066 {
00067 eps[i]=Mm->ip[ipp].strain[i];
00068 dpmax[i]=Mm->ip[ipp].eqother[i];
00069 }
00070
00071 if (fabs(eps[0])<s)
00072 sig[0]=taumax*pow(fabs(eps[0])/s,alfa)*sgn(eps[0]);
00073
00074 if ((fabs(eps[0])>s) && (fabs(eps[0])<ss))
00075 sig[0]=taumax*sgn(eps[0]);
00076
00077 if ((fabs(eps[0])>ss) && (fabs(eps[0])<sss))
00078 sig[0]=taumax-(taumax-tauf)/(sss-ss)*(eps[0]-ss);
00079
00080 if (fabs(eps[0])>sss)
00081 sig[0]=tauf*sgn(eps[0]);
00082
00083
00084 sig[1]=eps[1]*normal;
00085 if (sig[1]>fct)
00086 sig[1]=fct;
00087 if (sig[1]<-fcc)
00088 sig[1]=-fcc;
00089
00090 for (i=0;i<n;i++)
00091 Mm->ip[ipp].stress[i]=sig[i];
00092
00093 }
00094
00095 void cebfipcontactmat::updateval (long ipp, long im, long ido)
00096 {
00097 long i,n = Mm->givencompeqother(ipp, im);
00098
00099 for (i=0;i<n;i++){
00100 Mm->ip[ipp].eqother[ido+i]=Mm->ip[ipp].other[ido+i];
00101 }
00102 }
00103
00104 long cebfipcontactmat::sgn (double a)
00105 {
00106 long b;
00107 if (a==0) b=0;
00108 if (a<0) b=-1;
00109 if (a>0) b=1;
00110 return (b);
00111 }
00112