00001 #include <math.h>
00002 #include "chen.h"
00003 #include "global.h"
00004 #include "genfile.h"
00005 #include "intpoints.h"
00006 #include "vecttens.h"
00007 #include "alias.h"
00008 #include "stochdriver.h"
00009 #include "matrix.h"
00010
00011 chen::chen (void)
00012 {
00013 fyc=0.0;
00014 fyt=0.0;
00015 fybc=0.0;
00016 fc=0.0;
00017 ft=0.0;
00018 fbc=0.0;
00019 epsuc=0.0;
00020 epsut=0.0;
00021
00022 ay=0.0; au=0.0;
00023 ky=0.0; ku=0.0;
00024
00025
00026
00027
00028 state=2;
00029
00030
00031 }
00032
00033 chen::~chen (void)
00034 {
00035
00036 }
00037
00038
00039
00040
00041
00042
00043 void chen::read (FILE *in)
00044 {
00045 fscanf (in,"%lf %lf %lf %lf %lf %lf %lf %lf",&fyc,&fyt,&fybc,&fc,&ft,&fbc,&epsuc,&epsut);
00046 sra.read (in);
00047 }
00048
00049
00050
00051
00052
00053
00054
00055
00056
00057 void chen::matstiff (matrix &d,long ipp,long ido)
00058 {
00059 if (Mp->nlman->stmat==initial_stiff){
00060
00061 Mm->elmatstiff (d,ipp);
00062 }
00063 if (Mp->nlman->stmat==tangent_stiff){
00064
00065
00066
00067 matrix ad(d.m,d.n);
00068 Mm->elmatstiff (ad,ipp);
00069 tangentstiff (ad,d,ipp,ido);
00070 }
00071 }
00072
00073
00074 void chen::tangentstiff (matrix &d,matrix &td,long ipp,long ido)
00075 {
00076 long ncomp=Mm->ip[ipp].ncompstr;
00077 double denom1,denom2,denom,gamma;
00078 vector str,av(d.m),q(1),dfdq(1);
00079 matrix sig(3,3),am(d.m,d.n),dfds(3,3);
00080
00081 gamma=Mm->ip[ipp].eqother[ido+ncomp];
00082 if (gamma<1.0e-10){
00083 copym (d,td);
00084 }
00085 else{
00086 allocv (ncomp,str);
00087
00088
00089 Mm->givestress (0,ipp,str);
00090 vector_tensor (str,sig,Mm->ip[ipp].ssst,stress);
00091
00092 q[0]=Mm->ip[ipp].eqother[ido+ncomp+1];
00093
00094 deryieldfdsigma (sig,q,dfds);
00095 tensor_vector (str,dfds,Mm->ip[ipp].ssst,stress);
00096
00097 if (Mm->ip[ipp].ssst==planestress){
00098 vector auxstr(3);
00099 auxstr[0]=str[0];auxstr[1]=str[1];auxstr[2]=str[2];
00100 destrv (str);
00101 allocv (d.m,str);
00102 str[0]=auxstr[0];str[1]=auxstr[1];str[2]=auxstr[2];
00103 }
00104
00105 mxv (d,str,av);
00106 scprd (av,str,denom1);
00107
00108 scprd (str,str,denom2);
00109
00110 deryieldfdq (sig,q,dfdq);
00111
00112 denom2 = sqrt(denom2)*dfdq[0];
00113
00114 denom=denom1-denom2;
00115
00116 if (fabs(denom)<1.0e-10){
00117 copym (d,td);
00118 }
00119 else{
00120 vxv (str,str,am);
00121 mxm (d,am,td);
00122 mxm (td,d,am);
00123
00124 cmulm (1.0/denom,am);
00125
00126 subm (d,am,td);
00127 }
00128
00129 destrv (str);
00130 }
00131
00132 }
00133
00134
00135
00136
00137
00138
00139
00140
00141
00142 double chen::yieldfunction (matrix &sig,vector &q)
00143 {
00144 double f,invar,j2,zone,kappa;
00145 matrix dev(3,3);
00146
00147 deviator (sig,dev);
00148 invar=first_invar (sig);
00149 j2=second_invar(dev);
00150
00151 zone=sqrt(j2)+invar/sqrt(3.0);
00152
00153 if (state==1){
00154
00155 if(invar<0.0 && zone<0.0){
00156 ay = (fybc*fybc-fyc*fyc)/(2.0*fybc-fyc);
00157 au = (fbc*fbc-fc*fc)/(2.0*fbc-fc);
00158 ky = fyc*fybc*(2*fyc-fybc)/3.0/(2*fybc-fyc);
00159 ku = fc*fbc*(2*fc-fbc)/3.0/(2*fbc-fc);
00160
00161 if (q[0]>epsuc*epsuc) q[0]=epsuc*epsuc;
00162
00163 f = j2 + ay/3.0*invar - ky + q[0]*((au-ay)/epsuc/epsuc*invar/3.0 - (ku-ky)/epsuc/epsuc);
00164 }
00165
00166 else{
00167 ay = (fyc-fyt)/2.0;
00168 au = (fc-ft)/2.0;
00169 ky = fyc*fyt/6.0;
00170 ku = fc*ft/6.0;
00171
00172 if (q[0]>epsut*epsut) q[0]=epsut*epsut;
00173
00174 f = j2 - (invar*invar)/6.0 + ay/3.0*invar - ky +q[0]*((au-ay)/epsut/epsut*invar/3.0 - (ku-ky)/epsut/epsut);
00175 }
00176 }
00177 if (state==2){
00178 ay = (fybc*fybc-fyc*fyc)/(2.0*fybc-fyc);
00179 au = (fbc*fbc-fc*fc)/(2.0*fbc-fc);
00180 ky = fyc*fybc*(2*fyc-fybc)/3.0/(2*fybc-fyc);
00181 ku = fc*fbc*(2*fc-fbc)/3.0/(2*fbc-fc);
00182
00183 if (q[0]>epsuc*epsuc) q[0]=epsuc*epsuc;
00184
00185 f = j2 + ay/3.0*invar - ky + q[0]*((au-ay)/epsuc/epsuc*invar/3.0 - (ku-ky)/epsuc/epsuc);
00186 }
00187 if (state==3){
00188 ay = (fyc-fyt)/2.0;
00189 au = (fc-ft)/2.0;
00190 ky = fyc*fyt/6.0;
00191 ku = fc*ft/6.0;
00192
00193 if (q[0]>epsut*epsut) q[0]=epsut*epsut;
00194
00195 f = j2 - (invar*invar)/6.0 + ay/3.0*invar - ky +q[0]*((au-ay)/epsut/epsut*invar/3.0 - (ku-ky)/epsut/epsut);
00196 }
00197
00198 return f;
00199 }
00200
00201
00202
00203
00204
00205
00206
00207
00208 void chen::deryieldfdsigma (matrix &sig,vector &q,matrix &dfds)
00209 {
00210 double invar,j2,zone,kappa;
00211 matrix dev(3,3);
00212
00213 deviator (sig,dev);
00214 invar=first_invar (sig);
00215 j2=second_invar(dev);
00216
00217 zone=sqrt(j2)+invar/sqrt(3.0);
00218
00219 if (state==1){
00220
00221 if(invar<0.0 && zone<0.0){
00222 ay = (fybc*fybc-fyc*fyc)/(2.0*fybc-fyc);
00223 au = (fbc*fbc-fc*fc)/(2.0*fbc-fc);
00224 ky = fyc*fybc*(2*fyc-fybc)/3.0/(2*fybc-fyc);
00225 ku = fc*fbc*(2*fc-fbc)/3.0/(2*fbc-fc);
00226
00227 if (q[0]>epsuc*epsuc) q[0]=epsuc*epsuc;
00228
00229 dfds[0][0]=1.0/3.0*(2.0*sig[0][0]-sig[1][1]-sig[2][2]) + ay/3.0 + q[0]/3.0*(au-ay)/epsuc/epsuc;
00230 dfds[1][1]=1.0/3.0*(2.0*sig[1][1]-sig[0][0]-sig[2][2]) + ay/3.0 + q[0]/3.0*(au-ay)/epsuc/epsuc;
00231 dfds[2][2]=1.0/3.0*(2.0*sig[2][2]-sig[0][0]-sig[1][1]) + ay/3.0 + q[0]/3.0*(au-ay)/epsuc/epsuc;
00232
00233 dfds[0][1]=2*sig[0][1];
00234 dfds[1][0]=2*sig[0][1];
00235 dfds[0][2]=2*sig[0][2];
00236 dfds[2][0]=2*sig[0][2];
00237 dfds[1][2]=2*sig[1][2];
00238 dfds[2][1]=2*sig[1][2];
00239 }
00240
00241 else{
00242 ay = (fyc-fyt)/2.0;
00243 au = (fc-ft)/2.0;
00244 ky = fyc*fyt/6.0;
00245 ku = fc*ft/6.0;
00246
00247 if (q[0]>epsut*epsut) q[0]=epsut*epsut;
00248
00249 dfds[0][0]=1.0/3.0*(2.0*sig[0][0]-sig[1][1]-sig[2][2]) + ay/3.0 - 1.0/3.0*invar + q[0]/3.0*(au-ay)/epsut/epsut;
00250 dfds[1][1]=1.0/3.0*(2.0*sig[1][1]-sig[0][0]-sig[2][2]) + ay/3.0 - 1.0/3.0*invar + q[0]/3.0*(au-ay)/epsut/epsut;
00251 dfds[2][2]=1.0/3.0*(2.0*sig[2][2]-sig[0][0]-sig[1][1]) + ay/3.0 - 1.0/3.0*invar + q[0]/3.0*(au-ay)/epsut/epsut;
00252
00253 dfds[0][1]=2*sig[0][1];
00254 dfds[1][0]=2*sig[0][1];
00255 dfds[0][2]=2*sig[0][2];
00256 dfds[2][0]=2*sig[0][2];
00257 dfds[1][2]=2*sig[1][2];
00258 dfds[2][1]=2*sig[1][2];
00259 }
00260 }
00261 if (state==2){
00262 ay = (fybc*fybc-fyc*fyc)/(2.0*fybc-fyc);
00263 au = (fbc*fbc-fc*fc)/(2.0*fbc-fc);
00264 ky = fyc*fybc*(2*fyc-fybc)/3.0/(2*fybc-fyc);
00265 ku = fc*fbc*(2*fc-fbc)/3.0/(2*fbc-fc);
00266
00267 if (q[0]>epsuc*epsuc) q[0]=epsuc*epsuc;
00268
00269 dfds[0][0]=1.0/3.0*(2.0*sig[0][0]-sig[1][1]-sig[2][2]) + ay/3.0 + q[0]/3.0*(au-ay)/epsuc/epsuc;
00270 dfds[1][1]=1.0/3.0*(2.0*sig[1][1]-sig[0][0]-sig[2][2]) + ay/3.0 + q[0]/3.0*(au-ay)/epsuc/epsuc;
00271 dfds[2][2]=1.0/3.0*(2.0*sig[2][2]-sig[0][0]-sig[1][1]) + ay/3.0 + q[0]/3.0*(au-ay)/epsuc/epsuc;
00272
00273 dfds[0][1]=2*sig[0][1];
00274 dfds[1][0]=2*sig[0][1];
00275 dfds[0][2]=2*sig[0][2];
00276 dfds[2][0]=2*sig[0][2];
00277 dfds[1][2]=2*sig[1][2];
00278 dfds[2][1]=2*sig[1][2];
00279 }
00280 if (state==3){
00281 ay = (fyc-fyt)/2.0;
00282 au = (fc-ft)/2.0;
00283 ky = fyc*fyt/6.0;
00284 ku = fc*ft/6.0;
00285
00286 if (q[0]>epsut*epsut) q[0]=epsut*epsut;
00287
00288 dfds[0][0]=1.0/3.0*(2.0*sig[0][0]-sig[1][1]-sig[2][2]) + ay/3.0 - 1.0/3.0*invar + q[0]/3.0*(au-ay)/epsut/epsut;
00289 dfds[1][1]=1.0/3.0*(2.0*sig[1][1]-sig[0][0]-sig[2][2]) + ay/3.0 - 1.0/3.0*invar + q[0]/3.0*(au-ay)/epsut/epsut;
00290 dfds[2][2]=1.0/3.0*(2.0*sig[2][2]-sig[0][0]-sig[1][1]) + ay/3.0 - 1.0/3.0*invar + q[0]/3.0*(au-ay)/epsut/epsut;
00291
00292 dfds[0][1]=2*sig[0][1];
00293 dfds[1][0]=2*sig[0][1];
00294 dfds[0][2]=2*sig[0][2];
00295 dfds[2][0]=2*sig[0][2];
00296 dfds[1][2]=2*sig[1][2];
00297 dfds[2][1]=2*sig[1][2];
00298 }
00299 }
00300
00301
00302
00303
00304
00305
00306
00307 void chen::deryieldfdsigmadsigma (matrix &sig,matrix &dfdsds)
00308 {
00309 double invar,j2,zone;
00310 matrix dev(3,3);
00311
00312 deviator (sig,dev);
00313 invar=first_invar (sig);
00314 j2=second_invar(dev);
00315
00316 fillm (0.0,dfdsds);
00317
00318 zone = sqrt(j2)+invar/sqrt(3.0);
00319
00320 if (state==1){
00321
00322 if(invar<0.0 && zone<0.0){
00323 dfdsds[0][0] = 2.0/3.0;
00324 dfdsds[0][1] = -1.0/3.0;
00325 dfdsds[0][2] = -1.0/3.0;
00326
00327 dfdsds[1][0] = -1.0/3.0;
00328 dfdsds[1][1] = 2.0/3.0;
00329 dfdsds[1][2] = -1.0/3.0;
00330
00331 dfdsds[2][0] = -1.0/3.0;
00332 dfdsds[2][1] = -1.0/3.0;
00333 dfdsds[2][2] = 2.0/3.0;
00334
00335 dfdsds[3][3] = 2.0;
00336 dfdsds[4][4] = 2.0;
00337 dfdsds[5][5] = 2.0;
00338 }
00339
00340 else{
00341 dfdsds[0][0] = 1.0/3.0;
00342 dfdsds[0][1] = -2.0/3.0;
00343 dfdsds[0][2] = -2.0/3.0;
00344
00345 dfdsds[1][0] = -2.0/3.0;
00346 dfdsds[1][1] = 1.0/3.0;
00347 dfdsds[1][2] = -2.0/3.0;
00348
00349 dfdsds[2][0] = -2.0/3.0;
00350 dfdsds[2][1] = -2.0/3.0;
00351 dfdsds[2][2] = 1.0/3.0;
00352
00353 dfdsds[3][3] = 2.0;
00354 dfdsds[4][4] = 2.0;
00355 dfdsds[5][5] = 2.0;
00356 }
00357 }
00358 if (state==2){
00359 dfdsds[0][0] = 2.0/3.0;
00360 dfdsds[0][1] = -1.0/3.0;
00361 dfdsds[0][2] = -1.0/3.0;
00362
00363 dfdsds[1][0] = -1.0/3.0;
00364 dfdsds[1][1] = 2.0/3.0;
00365 dfdsds[1][2] = -1.0/3.0;
00366
00367 dfdsds[2][0] = -1.0/3.0;
00368 dfdsds[2][1] = -1.0/3.0;
00369 dfdsds[2][2] = 2.0/3.0;
00370
00371 dfdsds[3][3] = 2.0;
00372 dfdsds[4][4] = 2.0;
00373 dfdsds[5][5] = 2.0;
00374
00375 }
00376 if (state==3){
00377 dfdsds[0][0] = 1.0/3.0;
00378 dfdsds[0][1] = -2.0/3.0;
00379 dfdsds[0][2] = -2.0/3.0;
00380
00381 dfdsds[1][0] = -2.0/3.0;
00382 dfdsds[1][1] = 1.0/3.0;
00383 dfdsds[1][2] = -2.0/3.0;
00384
00385 dfdsds[2][0] = -2.0/3.0;
00386 dfdsds[2][1] = -2.0/3.0;
00387 dfdsds[2][2] = 1.0/3.0;
00388
00389 dfdsds[3][3] = 2.0;
00390 dfdsds[4][4] = 2.0;
00391 dfdsds[5][5] = 2.0;
00392 }
00393 }
00394
00395
00396
00397
00398
00399
00400
00401 void chen::deryieldfdq (matrix &sig,vector &q,vector &dfdq)
00402 {
00403 double invar,j2,zone;
00404 matrix dev(3,3);
00405
00406 deviator (sig,dev);
00407 invar=first_invar (sig);
00408 j2=second_invar(dev);
00409
00410 zone=sqrt(j2)+invar/sqrt(3.0);
00411
00412 if (state==1){
00413
00414 if(invar<0.0 && zone<0.0){
00415 ay = (fybc*fybc-fyc*fyc)/(2.0*fybc-fyc);
00416 au = (fbc*fbc-fc*fc)/(2.0*fbc-fc);
00417 ky = fyc*fybc*(2*fyc-fybc)/3.0/(2*fybc-fyc);
00418 ku = fc*fbc*(2*fc-fbc)/3.0/(2*fbc-fc);
00419
00420 dfdq[0] = (au-ay)/epsuc/epsuc*invar/3.0 - (ku-ky)/epsuc/epsuc;
00421 }
00422
00423 else{
00424 ay = (fyc-fyt)/2.0;
00425 au = (fc-ft)/2.0;
00426 ky = fyc*fyt/6.0;
00427 ku = fc*ft/6.0;
00428
00429 dfdq[0] = (au-ay)/epsut/epsut*invar/3.0 - (ku-ky)/epsut/epsut;
00430 }
00431 }
00432 if (state==2){
00433 ay = (fybc*fybc-fyc*fyc)/(2.0*fybc-fyc);
00434 au = (fbc*fbc-fc*fc)/(2.0*fbc-fc);
00435 ky = fyc*fybc*(2*fyc-fybc)/3.0/(2*fybc-fyc);
00436 ku = fc*fbc*(2*fc-fbc)/3.0/(2*fbc-fc);
00437
00438 dfdq[0] = (au-ay)/epsuc/epsuc*invar/3.0 - (ku-ky)/epsuc/epsuc;
00439 }
00440 if (state==3){
00441 ay = (fyc-fyt)/2.0;
00442 au = (fc-ft)/2.0;
00443 ky = fyc*fyt/6.0;
00444 ku = fc*ft/6.0;
00445
00446 dfdq[0] = (au-ay)/epsut/epsut*invar/3.0 - (ku-ky)/epsut/epsut;
00447 }
00448 }
00449
00450
00451
00452
00453 void chen::deryieldfdsigmadq (matrix &sig,vector &q,matrix &dfdsdq)
00454 {
00455 double invar,j2,zone;
00456 matrix dev(3,3);
00457
00458 deviator (sig,dev);
00459 invar=first_invar (sig);
00460 j2=second_invar(dev);
00461
00462 zone=sqrt(j2)+invar/sqrt(3.0);
00463
00464 if (state==1){
00465
00466 if(invar<0.0 && zone<0.0){
00467 ay = (fybc*fybc-fyc*fyc)/(2.0*fybc-fyc);
00468 au = (fbc*fbc-fc*fc)/(2.0*fbc-fc);
00469 ky = fyc*fybc*(2*fyc-fybc)/3.0/(2*fybc-fyc);
00470 ku = fc*fbc*(2*fc-fbc)/3.0/(2*fbc-fc);
00471
00472 dfdsdq[0][0]=(au-ay)/epsuc/epsuc/3.0;
00473 dfdsdq[1][0]=(au-ay)/epsuc/epsuc/3.0;
00474 dfdsdq[2][0]=(au-ay)/epsuc/epsuc/3.0;
00475
00476 dfdsdq[3][0]=0.0;
00477 dfdsdq[4][0]=0.0;
00478 dfdsdq[5][0]=0.0;
00479 }
00480
00481 else{
00482 ay = (fyc-fyt)/2.0;
00483 au = (fc-ft)/2.0;
00484 ky = fyc*fyt/6.0;
00485 ku = fc*ft/6.0;
00486
00487 dfdsdq[0][0]=(au-ay)/epsuc/epsuc/3.0;
00488 dfdsdq[1][0]=(au-ay)/epsuc/epsuc/3.0;
00489 dfdsdq[2][0]=(au-ay)/epsuc/epsuc/3.0;
00490
00491 dfdsdq[3][0]=0.0;
00492 dfdsdq[4][0]=0.0;
00493 dfdsdq[5][0]=0.0;
00494 }
00495 }
00496 if (state==2){
00497 ay = (fybc*fybc-fyc*fyc)/(2.0*fybc-fyc);
00498 au = (fbc*fbc-fc*fc)/(2.0*fbc-fc);
00499 ky = fyc*fybc*(2*fyc-fybc)/3.0/(2*fybc-fyc);
00500 ku = fc*fbc*(2*fc-fbc)/3.0/(2*fbc-fc);
00501
00502 dfdsdq[0][0]=(au-ay)/epsuc/epsuc/3.0;
00503 dfdsdq[1][0]=(au-ay)/epsuc/epsuc/3.0;
00504 dfdsdq[2][0]=(au-ay)/epsuc/epsuc/3.0;
00505
00506 dfdsdq[3][0]=0.0;
00507 dfdsdq[4][0]=0.0;
00508 dfdsdq[5][0]=0.0;
00509 }
00510 if (state==3){
00511 ay = (fyc-fyt)/2.0;
00512 au = (fc-ft)/2.0;
00513 ky = fyc*fyt/6.0;
00514 ku = fc*ft/6.0;
00515
00516 dfdsdq[0][0]=(au-ay)/epsuc/epsuc/3.0;
00517 dfdsdq[1][0]=(au-ay)/epsuc/epsuc/3.0;
00518 dfdsdq[2][0]=(au-ay)/epsuc/epsuc/3.0;
00519
00520 dfdsdq[3][0]=0.0;
00521 dfdsdq[4][0]=0.0;
00522 dfdsdq[5][0]=0.0;
00523 }
00524
00525 }
00526
00527 void chen::updateq(long ipp,double dgamma,vector &epsp,matrix &sig,vector &q)
00528 {
00529
00530 }
00531
00532
00533
00534
00535
00536
00537
00538
00539
00540
00541
00542 void chen::nlstresses (long ipp, long im, long ido)
00543 {
00544 long i,ni,n,nhard;
00545 double gamma,err;
00546
00547
00548 n = Mm->ip[ipp].ncompstr;
00549
00550 vector epsn(n),epsp(n),q(1);
00551
00552
00553 for (i=0;i<n;i++){
00554
00555 epsn[i]=Mm->ip[ipp].strain[i];
00556
00557 epsp[i]=Mm->ip[ipp].eqother[ido+i];
00558 }
00559
00560 gamma = Mm->ip[ipp].eqother[ido+n];
00561
00562 q[0]=Mm->ip[ipp].eqother[ido+n+1];
00563
00564
00565 if (sra.give_tsra () == cp){
00566 ni=sra.give_ni ();
00567 err=sra.give_err ();
00568
00569
00570
00571
00572 Mm->newton_stress_return_3 (ipp,im,ido,gamma,epsn,epsp,q,ni,err);
00573
00574 }
00575 else{
00576 fprintf (stderr,"\n\n wrong type of stress return algorithm is required in nlstresses (file %s, line %d).\n",__FILE__,__LINE__);
00577 abort ();
00578 }
00579
00580
00581
00582
00583 for (i=0;i<n;i++){
00584 Mm->ip[ipp].other[ido+i]=epsp[i];
00585 }
00586 Mm->ip[ipp].other[ido+n]=gamma;
00587 Mm->ip[ipp].other[ido+n+1]=q[0];
00588
00589 }
00590
00591 void chen::nonloc_nlstresses (long ipp, long im, long ido)
00592 {
00593 long i,ni, n = Mm->ip[ipp].ncompstr;
00594 double gamma,err;
00595 vector epsn(n),epsp(n),q(1);
00596
00597
00598 for (i=0;i<n;i++){
00599 epsn[i]=Mm->ip[ipp].strain[i];
00600 epsp[i]=Mm->ip[ipp].nonloc[i];
00601 }
00602 gamma = Mm->ip[ipp].eqother[ido+n];
00603 q[0] = Mm->ip[ipp].eqother[ido+n+1];
00604
00605
00606 if (sra.give_tsra () == cp){
00607 ni=sra.give_ni ();
00608 err=sra.give_err ();
00609
00610 Mm->cutting_plane (ipp,im,ido,gamma,epsn,epsp,q,ni,err);
00611 }
00612 else{
00613 fprintf (stderr,"\n\n wrong type of stress return algorithm is required in nlstresses (file %s, line %d).\n",__FILE__,__LINE__);
00614 abort ();
00615 }
00616
00617
00618 for (i=0;i<n;i++){
00619 Mm->ip[ipp].other[ido+i]=epsp[i];
00620 }
00621 Mm->ip[ipp].other[ido+n]=gamma;
00622 Mm->ip[ipp].other[ido+n+1]=q[0];
00623 }
00624
00625
00626
00627
00628
00629
00630 void chen::updateval (long ipp, long im, long ido)
00631 {
00632 long i,n = Mm->givencompeqother(ipp, im);
00633
00634 for (i=0;i<n;i++){
00635 Mm->ip[ipp].eqother[ido+i]=Mm->ip[ipp].other[ido+i];
00636 }
00637 }
00638
00639
00640
00641
00642
00643
00644 void chen::changeparam (atsel &atm,vector &val)
00645 {
00646 long i;
00647
00648 for (i=0;i<atm.num;i++){
00649 switch (atm.atrib[i]){
00650 case 0:{
00651 fyc=val[i];
00652 break;
00653 }
00654 case 1:{
00655 fyt=val[i];
00656 break;
00657 }
00658 case 2:{
00659 fybc=val[i];
00660 break;
00661 }
00662 case 3:{
00663 fc=val[i];
00664 break;
00665 }
00666 case 4:{
00667 ft=val[i];
00668 break;
00669 }
00670 case 5:{
00671 fbc=val[i];
00672 break;
00673 }
00674 default:{
00675 fprintf (stderr,"\n\n wrong number of atribute in function changeparam (file %s, line %d).\n",__FILE__,__LINE__);
00676 }
00677 }
00678 }
00679 }
00680
00681 double chen::give_consparam (long ipp,long ido)
00682 {
00683 long ncompstr;
00684 double gamma;
00685
00686 ncompstr=Mm->ip[ipp].ncompstr;
00687 gamma = Mm->ip[ipp].eqother[ido+ncompstr];
00688
00689 return gamma;
00690 }
00691
00692 long chen::give_num_interparam ()
00693 {
00694 return 1;
00695 }
00696
00697 void chen::give_interparam (long ipp,long ido,vector &q)
00698 {
00699 long ncompstr=Mm->ip[ipp].ncompstr;
00700
00701 q[0]=Mm->ip[ipp].eqother[ido+ncompstr+1];
00702 }