Changeset 21350


Ignore:
Timestamp:
11/07/16 17:22:11 (8 years ago)
Author:
adhikari
Message:

CHG: constants appearing in SLR rotational feedback double checked

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

Legend:

Unmodified
Added
Removed
  • TabularUnified issm/trunk-jpl/src/c/classes/Elements/Tria.cpp

    r21344 r21350  
    39303930                S=0; for(int i=0;i<NUMVERTICES;i++) S+=Sg_old[this->vertices[i]->Sid()]/NUMVERTICES;
    39313931               
     3932                /* Perturbation terms for moment of inertia (moi_list):
     3933                 * computed analytically (see Wu & Peltier, eqs 10 & 32)
     3934                 * also consistent with my GMD formulation!
     3935                 * ALL in geographic coordinates
     3936                 * */
    39323937                dI_list[0] = -4*PI*(rho_water*S*area)*pow(re,4)*(sin(late)*cos(late)*cos(longe))/eartharea;
    39333938                dI_list[1] = -4*PI*(rho_water*S*area)*pow(re,4)*(sin(late)*cos(late)*sin(longe))/eartharea;
    3934                 dI_list[2] = +4*PI*(rho_water*S*area)*pow(re,4)*(1-pow(cos(late),2))/eartharea;
     3939                dI_list[2] = +4*PI*(rho_water*S*area)*pow(re,4)*(1-pow(sin(late),2))/eartharea;
    39353940        }
    39363941        else if(this->inputs->Max(MaskIceLevelsetEnum)<0){
     
    39473952                dI_list[0] = -4*PI*(rho_ice*I*area)*pow(re,4)*(sin(late)*cos(late)*cos(longe))/eartharea;
    39483953                dI_list[1] = -4*PI*(rho_ice*I*area)*pow(re,4)*(sin(late)*cos(late)*sin(longe))/eartharea;
    3949                 dI_list[2] = +4*PI*(rho_ice*I*area)*pow(re,4)*(1-pow(cos(late),2))/eartharea;
     3954                dI_list[2] = +4*PI*(rho_ice*I*area)*pow(re,4)*(1-pow(sin(late),2))/eartharea;
    39503955        }
    39513956       
  • TabularUnified issm/trunk-jpl/src/c/classes/FemModel.cpp

    r21344 r21350  
    25392539        m1 = 1/(1-tide_love_k/fluid_love) * (1+load_love_k2)/(moi_p-moi_e) * moi_list[0];
    25402540        m2 = 1/(1-tide_love_k/fluid_love) * (1+load_love_k2)/(moi_p-moi_e) * moi_list[1];
    2541         m3 = -(1+load_love_k2)/moi_p * moi_list[2];
    2542 
     2541        m3 = -(1+load_love_k2)/moi_p * moi_list[2];     // term associated with fluid number (3-order-of-magnitude smaller) is negelected 
     2542
     2543        /* Green's function (1+k_2-h_2/g): checked against Glenn Milne's thesis Chapter 3 (eqs: 3.3-4, 3.10-11)
     2544         * Perturbation terms for angular velocity vector (m1, m2, m3): checked against Mitrovica (2005 Appendix) & Jensen et al (2013 Appendix A3)
     2545         * Sea level rotational feedback: checked against GMD eqs 8-9 (only first order terms, i.e., degree 2 order 0 & 1 considered)
     2546         * all DONE in Geographic coordinates: theta \in [-90,90], lambda \in [-180 180]
     2547         */
    25432548        for(int i=0;i<vertices->Size();i++){
    25442549                //Vertex* vertex=(Vertex*)vertices->GetObjectByOffset(i);
     
    25492554                /*only first order terms are considered now: */
    25502555                value=((1.0+tide_love_k-tide_love_h)/9.81)*pow(omega*radius[i],2.0)*
    2551                                                 (-m3/6.0 - m3*cos(2.0*lati)/2.0 - sin(lati)*cos(lati)*(m1*cos(longi)+m2*sin(longi)));
     2556                                                (-m3/6.0 + 0.5*m3*cos(2.0*lati) - 0.5*sin(2.*lati)*(m1*cos(longi)+m2*sin(longi)));
    25522557       
    25532558                pSgo_rot->SetValue(vertex->Sid(),value,INS_VAL); //INS_VAL ensures that you don't add several times
Note: See TracChangeset for help on using the changeset viewer.