Changeset 16915


Ignore:
Timestamp:
11/24/13 15:42:34 (11 years ago)
Author:
Mathieu Morlighem
Message:

CHG: moved more stuff from Tria to Element.cpp

Location:
issm/trunk-jpl/src/c/classes/Elements
Files:
7 edited

Legend:

Unmodified
Added
Removed

  • issm/trunk-jpl/src/c/classes/Elements/Element.cpp

    r16914 r16915  
    108108        return this->matpar->TMeltingPoint(pressure);
    109109}/*}}}*/
     110
     111void Element::ViscosityFS(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input){/*{{{*/
     112
     113        /*Intermediaries*/
     114        IssmDouble viscosity;
     115        IssmDouble epsilon3d[6]; /* epsilon=[exx,eyy,ezz,exy,exz,eyz];*/
     116        IssmDouble epsilon2d[3]; /* epsilon=[exx,eyy,exy];            */
     117
     118        if(vz_input){
     119                /*3D*/
     120                this->StrainRateFS(&epsilon3d[0],xyz_list,gauss,vx_input,vy_input,vz_input);
     121                material->GetViscosity3dFS(&viscosity, &epsilon3d[0]);
     122        }
     123        else{
     124                /*2D*/
     125                this->StrainRateSSA(&epsilon2d[0],xyz_list,gauss,vx_input,vy_input);
     126                material->GetViscosity2dvertical(&viscosity,&epsilon2d[0]);
     127        }
     128
     129        /*Assign output pointer*/
     130        *pviscosity=viscosity;
     131}
     132/*}}}*/
     133void Element::ViscosityL1L2(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* surface_input){/*{{{*/
     134        /*Compute the L1L2 viscosity
     135         *
     136         *      1
     137         * mu = - A^-1 (sigma'_e)^(1-n)
     138         *      2
     139         *
     140         * sigma'_e^2 = |sigma'_//|^2 + |sigma'_perp|^2 (see Perego 2012 eq. 17,18)
     141         *
     142         * L1L2 assumptions:
     143         *
     144         * (1) |eps_b|_// = A (|sigma'_//|^2 + |sigma'_perp|^2)^((n-1)/2) |sigma'_//|
     145         * (2) |sigma'_perp|^2 = |rho g (s-z) grad(s)|^2
     146         *
     147         * Assuming that n = 3, we have a polynom of degree 3 to solve (the only unkown is X=|sigma'_//|)
     148         *
     149         * A X^3 + A |rho g (s-z) grad(s)|^2 X - |eps_b|_// = 0     */
     150
     151        int        i;
     152        IssmDouble z,s,viscosity,p,q,delta;
     153        IssmDouble tau_perp,tau_par,eps_b,A;
     154        IssmDouble epsilonvx[5]; /*exx eyy exy exz eyz*/
     155        IssmDouble epsilonvy[5]; /*exx eyy exy exz eyz*/
     156        IssmDouble epsilon[5];   /*exx eyy exy exz eyz*/
     157        IssmDouble slope[3];
     158
     159        /*Check that both inputs have been found*/
     160        if (!vx_input || !vy_input || !surface_input) _error_("Input missing");
     161
     162        /*Get tau_perp*/
     163        surface_input->GetInputValue(&s,gauss);
     164        surface_input->GetInputDerivativeValue(&slope[0],xyz_list,gauss);
     165        z=GetZcoord(gauss);
     166        tau_perp = matpar->GetRhoIce() * matpar->GetG() * fabs(s-z)*sqrt(slope[0]*slope[0]+slope[1]*slope[1]);
     167
     168        /* Get eps_b*/
     169        vx_input->GetVxStrainRate3dHO(epsilonvx,xyz_list,gauss);
     170        vy_input->GetVyStrainRate3dHO(epsilonvy,xyz_list,gauss);
     171        for(i=0;i<5;i++) epsilon[i]=epsilonvx[i]+epsilonvy[i];
     172        eps_b = sqrt(epsilon[0]*epsilon[0] + epsilon[1]*epsilon[1] + epsilon[0]*epsilon[1] + epsilon[2]*epsilon[2]);
     173        if(eps_b==0.){
     174                *pviscosity = 2.5e+17;
     175                return;
     176        }
     177
     178        /*Get A*/
     179        _assert_(material->GetN()==3.0);
     180        A=material->GetA();
     181
     182        /*Solve for tau_perp (http://fr.wikipedia.org/wiki/Méthode_de_Cardan)*/
     183        p     = tau_perp *tau_perp;
     184        q     = - eps_b/A;
     185        delta = q *q + p*p*p*4./27.;
     186        _assert_(delta>0);
     187        tau_par = pow(0.5*(-q+sqrt(delta)),1./3.) - pow(0.5*(q+sqrt(delta)),1./3.);
     188
     189        /*Viscosity*/
     190        viscosity = 1./(2.*A*(tau_par*tau_par + tau_perp*tau_perp));
     191        _assert_(!xIsNan(viscosity));
     192        _assert_(viscosity > 0.);
     193
     194        /*Assign output pointer*/
     195        *pviscosity = viscosity;
     196}/*}}}*/
     197void Element::ViscosityHO(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input){/*{{{*/
     198
     199        /*Intermediaries*/
     200        IssmDouble viscosity;
     201        IssmDouble epsilon[5];/* epsilon=[exx,eyy,exy,exz,eyz];*/
     202
     203        this->StrainRateHO(&epsilon[0],xyz_list,gauss,vx_input,vy_input);
     204        material->GetViscosity3d(&viscosity, &epsilon[0]);
     205
     206        /*Assign output pointer*/
     207        *pviscosity=viscosity;
     208}/*}}}*/
     209void Element::ViscosityHODerivativeEpsSquare(IssmDouble* pmu_prime,IssmDouble* epsilon){/*{{{*/
     210        this->material->GetViscosityDerivativeEpsSquare(pmu_prime,epsilon);
     211}/*}}}*/
     212void Element::StrainRateHO(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input){/*{{{*/
     213        /*Compute the 3d Blatter/HOStrain Rate (5 components):
     214         *
     215         * epsilon=[exx eyy exy exz eyz]
     216         *
     217         * with exz=1/2 du/dz
     218         *      eyz=1/2 dv/dz
     219         *
     220         * the contribution of vz is neglected
     221         */
     222
     223        int i;
     224        IssmDouble epsilonvx[5];
     225        IssmDouble epsilonvy[5];
     226
     227        /*Check that both inputs have been found*/
     228        if (!vx_input || !vy_input){
     229                _error_("Input missing. Here are the input pointers we have for vx: " << vx_input << ", vy: " << vy_input << "\n");
     230        }
     231
     232        /*Get strain rate assuming that epsilon has been allocated*/
     233        vx_input->GetVxStrainRate3dHO(epsilonvx,xyz_list,gauss);
     234        vy_input->GetVyStrainRate3dHO(epsilonvy,xyz_list,gauss);
     235
     236        /*Sum all contributions*/
     237        for(i=0;i<5;i++) epsilon[i]=epsilonvx[i]+epsilonvy[i];
     238}/*}}}*/
     239void Element::ViscosityFSDerivativeEpsSquare(IssmDouble* pmu_prime,IssmDouble* epsilon){/*{{{*/
     240        this->material->GetViscosityDerivativeEpsSquare(pmu_prime,epsilon);
     241}/*}}}*/
     242void Element::StrainRateFS(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input){/*{{{*/
     243        /*Compute the 3d Strain Rate (6 components):
     244         *
     245         * epsilon=[exx eyy ezz exy exz eyz]
     246         */
     247
     248        IssmDouble epsilonvx[6];
     249        IssmDouble epsilonvy[6];
     250        IssmDouble epsilonvz[6];
     251
     252        /*Check that both inputs have been found*/
     253        if (!vx_input || !vy_input || !vz_input){
     254                _error_("Input missing. Here are the input pointers we have for vx: " << vx_input << ", vy: " << vy_input << ", vz: " << vz_input << "\n");
     255        }
     256
     257        /*Get strain rate assuming that epsilon has been allocated*/
     258        vx_input->GetVxStrainRate3d(epsilonvx,xyz_list,gauss);
     259        vy_input->GetVyStrainRate3d(epsilonvy,xyz_list,gauss);
     260        vz_input->GetVzStrainRate3d(epsilonvz,xyz_list,gauss);
     261
     262        /*Sum all contributions*/
     263        for(int i=0;i<6;i++) epsilon[i]=epsilonvx[i]+epsilonvy[i]+epsilonvz[i];
     264}/*}}}*/
     265void Element::ViscositySSA(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input){/*{{{*/
     266
     267        /*Intermediaries*/
     268        IssmDouble viscosity;
     269        IssmDouble epsilon[3];/* epsilon=[exx,eyy,exy];    */
     270
     271        this->StrainRateSSA(&epsilon[0],xyz_list,gauss,vx_input,vy_input);
     272        material->GetViscosity2d(&viscosity, &epsilon[0]);
     273
     274        /*Assign output pointer*/
     275        *pviscosity=viscosity;
     276}/*}}}*/
     277void Element::ViscositySSADerivativeEpsSquare(IssmDouble* pmu_prime,IssmDouble* epsilon){/*{{{*/
     278        this->material->GetViscosity2dDerivativeEpsSquare(pmu_prime,epsilon);
     279}/*}}}*/
     280void Element::StrainRateSSA(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input){/*{{{*/
     281
     282        int i;
     283        IssmDouble epsilonvx[3];
     284        IssmDouble epsilonvy[3];
     285
     286        /*Check that both inputs have been found*/
     287        if (!vx_input || !vy_input){
     288                _error_("Input missing. Here are the input pointers we have for vx: " << vx_input << ", vy: " << vy_input << "\n");
     289        }
     290
     291        /*Get strain rate assuming that epsilon has been allocated*/
     292        vx_input->GetVxStrainRate2d(epsilonvx,xyz_list,gauss);
     293        vy_input->GetVyStrainRate2d(epsilonvy,xyz_list,gauss);
     294
     295        /*Sum all contributions*/
     296        for(i=0;i<3;i++) epsilon[i]=epsilonvx[i]+epsilonvy[i];
     297}/*}}}*/

  • issm/trunk-jpl/src/c/classes/Elements/Element.h

    r16914 r16915  
    6262                IssmDouble GetMaterialParameter(int enum_in);
    6363                bool       IsFloating();
     64                void       StrainRateSSA(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input);
     65                void       StrainRateHO(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input);
     66                void       StrainRateFS(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input);
    6467                IssmDouble TMeltingPoint(IssmDouble pressure);
     68                void       ViscosityFS(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input);
     69                void       ViscosityHO(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input);
     70                void       ViscosityL1L2(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* surf);
     71                void       ViscositySSA(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input);
     72                void       ViscositySSADerivativeEpsSquare(IssmDouble* pmu_prime,IssmDouble* epsilon);
     73                void       ViscosityHODerivativeEpsSquare(IssmDouble* pmu_prime,IssmDouble* epsilon);
     74                void       ViscosityFSDerivativeEpsSquare(IssmDouble* pmu_prime,IssmDouble* epsilon);
    6575
    6676                /*Virtual functions*/
    … …  
    186196                virtual void   ValueP1DerivativesOnGauss(IssmDouble* dvalue,IssmDouble* values,IssmDouble* xyz_list,Gauss* gauss)=0;
    187197                virtual void   ViscousHeating(IssmDouble* pphi,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input)=0;
    188                 virtual void   ViscosityFS(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input)=0;
    189                 virtual void   ViscosityHO(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input)=0;
    190                 virtual void   ViscosityL1L2(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* surf)=0;
    191                 virtual void   ViscositySSA(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input)=0;
    192                 virtual void   ViscositySSADerivativeEpsSquare(IssmDouble* pmu_prime,IssmDouble* epsilon)=0;
    193                 virtual void   ViscosityHODerivativeEpsSquare(IssmDouble* pmu_prime,IssmDouble* epsilon)=0;
    194                 virtual void   ViscosityFSDerivativeEpsSquare(IssmDouble* pmu_prime,IssmDouble* epsilon)=0;
    195                 virtual void   StrainRateSSA(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input)=0;
    196                 virtual void   StrainRateHO(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input)=0;
    197                 virtual void   StrainRateFS(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input)=0;
     198
    198199                virtual int    VelocityInterpolation()=0;
    199200                virtual int    PressureInterpolation()=0;

  • issm/trunk-jpl/src/c/classes/Elements/Penta.cpp

    r16914 r16915  
    34523452        /*Assign output pointer*/
    34533453        *pphi = phi;
    3454 }
    3455 /*}}}*/
    3456 /*FUNCTION Penta::ViscosityFS{{{*/
    3457 void Penta::ViscosityFS(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input){
    3458 
    3459         /*Intermediaries*/
    3460         IssmDouble viscosity;
    3461         IssmDouble epsilon[6]; /* epsilon=[exx,eyy,ezz,exy,exz,eyz];*/
    3462 
    3463         this->StrainRateFS(&epsilon[0],xyz_list,gauss,vx_input,vy_input,vz_input);
    3464         material->GetViscosity3dFS(&viscosity, &epsilon[0]);
    3465 
    3466         /*Assign output pointer*/
    3467         *pviscosity=viscosity;
    3468 }
    3469 /*}}}*/
    3470 /*FUNCTION Penta::ViscosityL1L2{{{*/
    3471 void Penta::ViscosityL1L2(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* surface_input){
    3472         /*Compute the L1L2 viscosity
    3473          *
    3474          *      1
    3475          * mu = - A^-1 (sigma'_e)^(1-n)
    3476          *      2
    3477          *
    3478          * sigma'_e^2 = |sigma'_//|^2 + |sigma'_perp|^2 (see Perego 2012 eq. 17,18)
    3479          *
    3480          * L1L2 assumptions:
    3481          *
    3482          * (1) |eps_b|_// = A (|sigma'_//|^2 + |sigma'_perp|^2)^((n-1)/2) |sigma'_//|
    3483          * (2) |sigma'_perp|^2 = |rho g (s-z) grad(s)|^2
    3484          *
    3485          * Assuming that n = 3, we have a polynom of degree 3 to solve (the only unkown is X=|sigma'_//|)
    3486          *
    3487          * A X^3 + A |rho g (s-z) grad(s)|^2 X - |eps_b|_// = 0     */
    3488 
    3489         int        i;
    3490         IssmDouble z,s,viscosity,p,q,delta;
    3491         IssmDouble tau_perp,tau_par,eps_b,A;
    3492         IssmDouble epsilonvx[5]; /*exx eyy exy exz eyz*/
    3493         IssmDouble epsilonvy[5]; /*exx eyy exy exz eyz*/
    3494         IssmDouble epsilon[5];   /*exx eyy exy exz eyz*/
    3495         IssmDouble z_list[NUMVERTICES];
    3496         IssmDouble slope[3];
    3497 
    3498         /*Check that both inputs have been found*/
    3499         if (!vx_input || !vy_input || !surface_input) _error_("Input missing");
    3500 
    3501         /*Get tau_perp*/
    3502         for(i=0;i<NUMVERTICES;i++) z_list[i]=xyz_list[3*i+2];
    3503         surface_input->GetInputValue(&s,gauss);
    3504         surface_input->GetInputDerivativeValue(&slope[0],xyz_list,gauss);
    3505         PentaRef::GetInputValue(&z,&z_list[0],gauss);
    3506         tau_perp = matpar->GetRhoIce() * matpar->GetG() * fabs(s-z)*sqrt(slope[0]*slope[0]+slope[1]*slope[1]);
    3507 
    3508         /* Get eps_b*/
    3509         vx_input->GetVxStrainRate3dHO(epsilonvx,xyz_list,gauss);
    3510         vy_input->GetVyStrainRate3dHO(epsilonvy,xyz_list,gauss);
    3511         for(i=0;i<5;i++) epsilon[i]=epsilonvx[i]+epsilonvy[i];
    3512         eps_b = sqrt(epsilon[0]*epsilon[0] + epsilon[1]*epsilon[1] + epsilon[0]*epsilon[1] + epsilon[2]*epsilon[2]);
    3513         if(eps_b==0.){
    3514                 *pviscosity = 2.5e+17;
    3515                 return;
    3516         }
    3517 
    3518         /*Get A*/
    3519         _assert_(material->GetN()==3.0);
    3520         A=material->GetA();
    3521 
    3522         /*Solve for tau_perp (http://fr.wikipedia.org/wiki/Méthode_de_Cardan)*/
    3523         p     = tau_perp *tau_perp;
    3524         q     = - eps_b/A;
    3525         delta = q *q + p*p*p*4./27.;
    3526         _assert_(delta>0);
    3527         tau_par = pow(0.5*(-q+sqrt(delta)),1./3.) - pow(0.5*(q+sqrt(delta)),1./3.);
    3528 
    3529         /*Viscosity*/
    3530         viscosity = 1./(2.*A*(tau_par*tau_par + tau_perp*tau_perp));
    3531         _assert_(!xIsNan(viscosity));
    3532         _assert_(viscosity > 0.);
    3533 
    3534         /*Assign output pointer*/
    3535         *pviscosity = viscosity;
    3536 }
    3537 /*}}}*/
    3538 /*FUNCTION Penta::ViscosityHO{{{*/
    3539 void Penta::ViscosityHO(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input){
    3540 
    3541         /*Intermediaries*/
    3542         IssmDouble viscosity;
    3543         IssmDouble epsilon[5];/* epsilon=[exx,eyy,exy,exz,eyz];*/
    3544 
    3545         this->StrainRateHO(&epsilon[0],xyz_list,gauss,vx_input,vy_input);
    3546         material->GetViscosity3d(&viscosity, &epsilon[0]);
    3547 
    3548         /*Assign output pointer*/
    3549         *pviscosity=viscosity;
    3550 }
    3551 /*}}}*/
    3552 /*FUNCTION Penta::ViscosityHODerivativeEpsSquare{{{*/
    3553 void Penta::ViscosityHODerivativeEpsSquare(IssmDouble* pmu_prime,IssmDouble* epsilon){
    3554         this->material->GetViscosityDerivativeEpsSquare(pmu_prime,epsilon);
    3555 }
    3556 /*}}}*/
    3557 /*FUNCTION Penta::StrainRateHO{{{*/
    3558 void Penta::StrainRateHO(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input){
    3559         /*Compute the 3d Blatter/HOStrain Rate (5 components):
    3560          *
    3561          * epsilon=[exx eyy exy exz eyz]
    3562          *
    3563          * with exz=1/2 du/dz
    3564          *      eyz=1/2 dv/dz
    3565          *
    3566          * the contribution of vz is neglected
    3567          */
    3568 
    3569         int i;
    3570         IssmDouble epsilonvx[5];
    3571         IssmDouble epsilonvy[5];
    3572 
    3573         /*Check that both inputs have been found*/
    3574         if (!vx_input || !vy_input){
    3575                 _error_("Input missing. Here are the input pointers we have for vx: " << vx_input << ", vy: " << vy_input << "\n");
    3576         }
    3577 
    3578         /*Get strain rate assuming that epsilon has been allocated*/
    3579         vx_input->GetVxStrainRate3dHO(epsilonvx,xyz_list,gauss);
    3580         vy_input->GetVyStrainRate3dHO(epsilonvy,xyz_list,gauss);
    3581 
    3582         /*Sum all contributions*/
    3583         for(i=0;i<5;i++) epsilon[i]=epsilonvx[i]+epsilonvy[i];
    3584 }
    3585 /*}}}*/
    3586 /*FUNCTION Penta::ViscosityFSDerivativeEpsSquare{{{*/
    3587 void Penta::ViscosityFSDerivativeEpsSquare(IssmDouble* pmu_prime,IssmDouble* epsilon){
    3588         this->material->GetViscosityDerivativeEpsSquare(pmu_prime,epsilon);
    3589 }
    3590 /*}}}*/
    3591 /*FUNCTION Penta::StrainRateFS{{{*/
    3592 void Penta::StrainRateFS(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input){
    3593         /*Compute the 3d Strain Rate (6 components):
    3594          *
    3595          * epsilon=[exx eyy ezz exy exz eyz]
    3596          */
    3597 
    3598         IssmDouble epsilonvx[6];
    3599         IssmDouble epsilonvy[6];
    3600         IssmDouble epsilonvz[6];
    3601 
    3602         /*Check that both inputs have been found*/
    3603         if (!vx_input || !vy_input || !vz_input){
    3604                 _error_("Input missing. Here are the input pointers we have for vx: " << vx_input << ", vy: " << vy_input << ", vz: " << vz_input << "\n");
    3605         }
    3606 
    3607         /*Get strain rate assuming that epsilon has been allocated*/
    3608         vx_input->GetVxStrainRate3d(epsilonvx,xyz_list,gauss);
    3609         vy_input->GetVyStrainRate3d(epsilonvy,xyz_list,gauss);
    3610         vz_input->GetVzStrainRate3d(epsilonvz,xyz_list,gauss);
    3611 
    3612         /*Sum all contributions*/
    3613         for(int i=0;i<6;i++) epsilon[i]=epsilonvx[i]+epsilonvy[i]+epsilonvz[i];
    36143454}
    36153455/*}}}*/

  • issm/trunk-jpl/src/c/classes/Elements/Penta.h

    r16914 r16915  
    278278                void           TransformStiffnessMatrixCoord(ElementMatrix* Ke,int numnodes,int* transformenum_list){_error_("not implemented yet");};
    279279                void           TransformStiffnessMatrixCoord(ElementMatrix* Ke,Node** nodes_list,int numnodes,int* transformenum_list);
    280                 void           ViscosityFS(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input);
    281                 void           ViscosityHO(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input);
    282                 void           ViscosityL1L2(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* surf);
    283                 void           ViscositySSA(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input){_error_("not implemented");};
    284                 void           ViscosityHODerivativeEpsSquare(IssmDouble* pmu_prime,IssmDouble* epsilon);
    285                 void           ViscosityFSDerivativeEpsSquare(IssmDouble* pmu_prime,IssmDouble* epsilon);
    286                 void           ViscositySSADerivativeEpsSquare(IssmDouble* pmu_prime,IssmDouble* epsilon){_error_("not implemented");};
    287                 void           StrainRateFS(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input);
    288                 void           StrainRateHO(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input);
    289                 void           StrainRateSSA(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input){_error_("not implemented");};
    290280
    291281                #ifdef _HAVE_STRESSBALANCE_

  • issm/trunk-jpl/src/c/classes/Elements/Seg.h

    r16914 r16915  
    152152                void        VerticalSegmentIndices(int** pindices,int* pnumseg){_error_("not implemented yet");};
    153153                void        ViscousHeating(IssmDouble* pphi,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input){_error_("not implemented yet");};
    154                 void        ViscosityFS(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input){_error_("not implemented");};
    155                 void        ViscosityHO(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input){_error_("not implemented");};
    156                 void        ViscosityL1L2(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* surf){_error_("not implemented yet");};
    157                 void        ViscositySSA(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input){_error_("not implemented");};
    158                 void        ViscosityHODerivativeEpsSquare(IssmDouble* pmu_prime,IssmDouble* epsilon){_error_("not implemented");};
    159                 void        ViscosityFSDerivativeEpsSquare(IssmDouble* pmu_prime,IssmDouble* epsilon){_error_("not implemented");};
    160                 void        ViscositySSADerivativeEpsSquare(IssmDouble* pmu_prime,IssmDouble* epsilon){_error_("not implemented");};
    161                 void        StrainRateFS(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input){_error_("not implemented");};
    162                 void        StrainRateHO(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input){_error_("not implemented");};
    163                 void        StrainRateSSA(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input){_error_("not implemented");};
    164154                bool        IsZeroLevelset(int levelset_enum){_error_("not implemented");};
    165155                void        ZeroLevelsetCoordinates(IssmDouble** pxyz_zero,IssmDouble* xyz_list,int levelsetenum){_error_("not implemented");};

  • issm/trunk-jpl/src/c/classes/Elements/Tria.cpp

    r16914 r16915  
    30003000        _assert_(this->vertices);
    30013001        return this->vertices[vertexindex]->Connectivity();
    3002 }
    3003 /*}}}*/
    3004 /*FUNCTION Tria::ViscosityFS{{{*/
    3005 void Tria::ViscosityFS(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input){
    3006 
    3007         /*Intermediaries*/
    3008         IssmDouble viscosity;
    3009         IssmDouble epsilon[3];/* epsilon=[exx,eyy,exy];    */
    3010 
    3011         this->StrainRateSSA(&epsilon[0],xyz_list,gauss,vx_input,vy_input);
    3012         material->GetViscosity2dvertical(&viscosity, &epsilon[0]);
    3013 
    3014         /*Assign output pointer*/
    3015         *pviscosity=viscosity;
    3016 }
    3017 /*}}}*/
    3018 /*FUNCTION Tria::ViscositySSA{{{*/
    3019 void Tria::ViscositySSA(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input){
    3020 
    3021         /*Intermediaries*/
    3022         IssmDouble viscosity;
    3023         IssmDouble epsilon[3];/* epsilon=[exx,eyy,exy];    */
    3024 
    3025         this->StrainRateSSA(&epsilon[0],xyz_list,gauss,vx_input,vy_input);
    3026         material->GetViscosity2d(&viscosity, &epsilon[0]);
    3027 
    3028         /*Assign output pointer*/
    3029         *pviscosity=viscosity;
    3030 }
    3031 /*}}}*/
    3032 /*FUNCTION Tria::ViscositySSADerivativeEpsSquare{{{*/
    3033 void Tria::ViscositySSADerivativeEpsSquare(IssmDouble* pmu_prime,IssmDouble* epsilon){
    3034         this->material->GetViscosity2dDerivativeEpsSquare(pmu_prime,epsilon);
    3035 }
    3036 /*}}}*/
    3037 /*FUNCTION Tria::StrainRateSSA{{{*/
    3038 void Tria::StrainRateSSA(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input){
    3039 
    3040         int i;
    3041         IssmDouble epsilonvx[3];
    3042         IssmDouble epsilonvy[3];
    3043 
    3044         /*Check that both inputs have been found*/
    3045         if (!vx_input || !vy_input){
    3046                 _error_("Input missing. Here are the input pointers we have for vx: " << vx_input << ", vy: " << vy_input << "\n");
    3047         }
    3048 
    3049         /*Get strain rate assuming that epsilon has been allocated*/
    3050         vx_input->GetVxStrainRate2d(epsilonvx,xyz_list,gauss);
    3051         vy_input->GetVyStrainRate2d(epsilonvy,xyz_list,gauss);
    3052 
    3053         /*Sum all contributions*/
    3054         for(i=0;i<3;i++) epsilon[i]=epsilonvx[i]+epsilonvy[i];
    30553002}
    30563003/*}}}*/

  • issm/trunk-jpl/src/c/classes/Elements/Tria.h

    r16914 r16915  
    309309                void           TransformStiffnessMatrixCoord(ElementMatrix* Ke,Node** nodes_list,int numnodes,int* transformenum_list){_error_("not implemented yet");};
    310310                void           ViscousHeating(IssmDouble* pphi,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input){_error_("not implemented yet");};
    311                 void           ViscosityFS(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input);
    312                 void           ViscosityHO(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input){_error_("not implemented yet");};
    313                 void           ViscosityL1L2(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* surf){_error_("not implemented yet");};
    314                 void           ViscositySSA(IssmDouble* pviscosity,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input);
    315                 void           ViscosityHODerivativeEpsSquare(IssmDouble* pmu_prime,IssmDouble* epsilon){_error_("not implemented");};
    316                 void           ViscosityFSDerivativeEpsSquare(IssmDouble* pmu_prime,IssmDouble* epsilon){_error_("not implemented");};
    317                 void           ViscositySSADerivativeEpsSquare(IssmDouble* pmu_prime,IssmDouble* epsilon);
    318                 void           StrainRateFS(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input,Input* vz_input){_error_("not implemented");};
    319                 void           StrainRateHO(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input){_error_("not implemented");};
    320                 void           StrainRateSSA(IssmDouble* epsilon,IssmDouble* xyz_list,Gauss* gauss,Input* vx_input,Input* vy_input);
    321311
    322312                #ifdef _HAVE_STRESSBALANCE_
Note: See TracChangeset for help on using the changeset viewer.