Changeset 2959 in Sophya for trunk

Timestamp:

Jun 1, 2006, 4:34:31 PM (19 years ago)

Author:

ansari

Message:

correction fonctions optimisees calcul Ylm en version <r_4> - Reza 1/06/2006

File:

• trunk/SophyaLib/Samba/lambdaBuilder.cc (modified) (2 diffs)

trunk/SophyaLib/Samba/lambdaBuilder.cc

@@ -223 +223 @@
       
       //DEL        lambda_(l,m) = lam_2*lam_mm;
-      *almp += (lam_2*lam_mm) * phi;
-      almp++;
+      *almp += (lam_2*lam_mm) * phi;    almp++;
       
       b_rec=1./a_rec;
@@ -244 +243 @@
                                         TVector< complex<r_4> >& phase, Alm<r_4> & alm) 
 {
-  TVector< complex<r_8> > phase8;
-  phase8 = phase;
-  Alm<r_8> alm8(alm.Lmax());
-  ComputeAlmFrPhase(theta, lmax, mmax, phase8, alm8);
-  for(sa_size_t k=0; k<alm.Size(); k++) 
-    alm(k) = complex<r_4>((r_4)alm8(k).real() , (r_4)alm8(k).imag());
-  return;
+  updateArrayRecurrence(lmax);
+  r_8 cth = cos(theta);
+  r_8 sth = sin(theta);
+  
+  r_8 bignorm2 = 1.e268; // = 1e-20*1.d288
+  
+  r_8 lam_mm = 1. / sqrt(4.*Pi) *bignorm2;
+  r_8 lam_0, lam_1, lam_2;
+  register r_8 a_rec, b_rec;
+  register r_4 cff;
+  complex<r_4>* almp = alm.columnData(0); 
+  r_8* arecp = a_recurrence_.columnData(0); 
+
+  int_4 l, m;
+  for (m=0; m<=mmax;m++)  {
+    //MOV^ complex<r_4>* almp = alm.columnData(m); 
+    complex<r_4>  phi = phase(m);
+    
+    cff = lam_mm / bignorm2;
+    *almp += cff*phi;  almp++;
+
+    lam_0=0.;
+    lam_1=1. /bignorm2 ;
+    
+    //MOV^ r_8* arecp = a_recurrence_.columnData(m); 
+    a_rec = *arecp; arecp++; // a_recurrence_(m,m);
+    b_rec = 0.;
+    for (l=m+1; l<=lmax; l++)   {
+      lam_2 = (cth*lam_1-b_rec*lam_0)*a_rec;
+      
+      //DEL        lambda_(l,m) = lam_2*lam_mm;
+      cff = (lam_2*lam_mm);
+      *almp += cff * phi;    almp++;
+      
+      b_rec=1./a_rec;
+      //           a_rec= LWK->a_recurr(l,m);
+      a_rec= *arecp; arecp++; // a_recurrence_(l,m);
+      lam_0 = lam_1;
+      lam_1 = lam_2;
+    }
+    lam_mm = -lam_mm*sth* sqrt( (2.*m+3.)/ (2.*m+2.) );
+  }
 }