Changeset 20054


Ignore:
Timestamp:
02/02/16 19:18:00 (9 years ago)
Author:
Eric.Larour
Message:

CHG: faster parallel way of precomputing legendre * love polynomial values.
Fixed huge bug for AD! Pn, Pn1 and Pn2 were not getting renewed!

File:
1 edited

Legend:

Unmodified
Added
Removed
  • issm/trunk-jpl/src/c/analyses/SealevelriseAnalysis.cpp

    r20044 r20054  
    7373               
    7474                /*compute combined legendre + love number (elastic green function:*/
    75                 #ifdef _HAVE_ADOLC_
    76                 for(int i=0;i<M;i++){
    77                         IssmDouble alpha,x;
    78                         alpha= reCast<IssmDouble>(i)*degacc * PI / 180.0;
    79 
    80                         G_elastic[i]= (love_k[nl-1]-love_h[nl-1])/2.0/sin(alpha/2.0);
    81                         for (int n=0;n<nl;n++) {
    82                                 IssmDouble Pn,Pn1,Pn2;
    83                                 IssmDouble deltalove;
    84 
    85                                 deltalove = (love_k[n]-love_k[nl-1]-love_h[n]+love_h[nl-1]);
    86 
    87                                 if(n==0)Pn=1;
    88                                 else if(n==1)Pn=cos(alpha);
    89                                 else Pn= ( (2*n-1)*cos(alpha)*Pn1 - (n-1)*Pn2 ) /n;
    90                                 Pn2=Pn1; Pn1=Pn;
    91 
    92                                 G_elastic[i] += deltalove*Pn;
    93                         }
    94                 }
    95                 #else
    9675                m=DetermineLocalSize(M,IssmComm::GetComm());
    9776                GetOwnershipBoundariesFromRange(&lower_row,&upper_row,m,IssmComm::GetComm());
     
    10382
    10483                        G_elastic_local[i-lower_row]= (love_k[nl-1]-love_h[nl-1])/2.0/sin(alpha/2.0);
     84                        IssmDouble Pn,Pn1,Pn2;
    10585                        for (int n=0;n<nl;n++) {
    106                                 IssmDouble Pn,Pn1,Pn2;
    10786                                IssmDouble deltalove;
    10887
     
    135114
    136115                /*}}}*/
    137                 #endif
    138116
    139117                /*Avoid singularity at 0: */
Note: See TracChangeset for help on using the changeset viewer.