Changeset 2423 in Sophya


Ignore:
Timestamp:
Aug 28, 2003, 3:29:03 PM (22 years ago)
Author:
cmv
Message:

Quelques ameliorations mineures cmv 28/8/2003

Location:
trunk/SophyaLib/NTools
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • trunk/SophyaLib/NTools/integ.cc

    r2283 r2423  
    117117
    118118Integrator::~Integrator()
    119 {delete mFunc;}
     119{if(mFunc) delete mFunc;}
    120120
    121121//++
     
    222222  Limits(xmin,xmax);
    223223  return Value();
     224}
     225
     226void
     227Integrator::Print(int lp)
     228{
     229 cout<<"Integrator between "<<mXMin<<" and "<<mXMax
     230     <<" in "<<mNStep<<" steps"<<endl;
     231 if(lp>1) cout<<"mFunc="<<mFunc<<"mGFF="<<mGFF
     232              <<"mReqPrec="<<mReqPrec<<"mDX="<<mDX<<endl;
    224233}
    225234
     
    332341//      fonction aux zéros des polynomes de Legendre. Avec les poids
    333342//      qui vont bien, GL d'ordre n est exacte pour des polynomes de
    334 //      degré <= 2n+1.
     343//      degré <= 2n+1 (monome le + haut x^(2*n-1).
    335344//      Impossible de demander une précision donnée.
    336345//
     
    408417GLInteg::~GLInteg()
    409418{
    410   delete[] mXPos;
    411   delete[] mWeights;
    412   mXPos = mWeights = NULL;
     419 if(mXPos)    delete[] mXPos;     mXPos    = NULL;
     420 if(mWeights) delete[] mWeights;  mWeights = NULL;
    413421}
    414422
     
    416424GLInteg::InvalWeights()
    417425{
    418   delete[] mXPos;
    419   delete[] mWeights;
    420   mXPos = mWeights = NULL;
     426 if(mXPos)    delete[] mXPos;     mXPos    = NULL;
     427 if(mWeights) delete[] mWeights;  mWeights = NULL;
    421428}
    422429
     
    425432GLInteg::LimitsChanged()
    426433{
    427   ComputeBounds();
     434 ComputeBounds();
    428435}
    429436
     
    431438GLInteg::StepsChanged()
    432439{
    433   ComputeBounds();
     440 ComputeBounds();
    434441}
    435442
     
    437444GLInteg::AddBound(double x)
    438445{
    439   if (x<=mXMin || x>=mXMax) THROW(rangeCheckErr);
    440   // On introduira les classes d'exections apres reflexion et de maniere systematique (Rz+cmv)
    441   // if (x<=mXMin || x>=mXMax) throw range_error("GLInteg::AddBound  bound outside interval");
    442   mBounds.insert(x);
    443   return *this;
     446 if (x<=mXMin || x>=mXMax) THROW(rangeCheckErr);
     447 // On introduira les classes d'exections apres reflexion et de maniere systematique (Rz+cmv)
     448 // if (x<=mXMin || x>=mXMax) throw range_error("GLInteg::AddBound  bound outside interval");
     449 mBounds.insert(x);
     450 return *this;
    444451}
    445452
     
    448455GLInteg::SetOrder(int order)
    449456{
    450   mOrder = order;
    451   InvalWeights();
    452   return *this;
     457 mOrder = order;
     458 InvalWeights();
     459 return *this;
    453460}
    454461
     
    458465GLInteg::Value()
    459466{
    460   if (!mXPos) ComputeWeights();
    461   double s = 0;
    462   set<double>::iterator i=mBounds.begin();
    463   set<double>::iterator j=i; j++;
    464   while (j != mBounds.end()) {
    465     double s1 = 0;
    466     double x1 = *i;
    467     double x2 = *j;
    468     for (int k=0; k<mOrder; k++)
    469        s1 += mWeights[k] * FVal(x1 + (x2-x1)*mXPos[k]);
    470     s += s1*(x2-x1);
    471     i++; j++;
    472   }
    473   return s;
     467 if (!mXPos) ComputeWeights();
     468 double s = 0;
     469 set<double>::iterator i=mBounds.begin();
     470 set<double>::iterator j=i; j++;
     471 while (j != mBounds.end()) {
     472   double s1 = 0;
     473   double x1 = *i;
     474   double x2 = *j;
     475   for(int k=0; k<mOrder; k++)
     476      s1 += mWeights[k] * FVal(x1 + (x2-x1)*mXPos[k]);
     477   s += s1*(x2-x1);
     478   i++; j++;
     479 }
     480 return s;
    474481}
    475482
     
    485492GLInteg::ComputeWeights()
    486493{
    487    const double EPS_gauleg = 3.0e-11;
     494   const double EPS_gauleg = 5.0e-12;
    488495
    489496   int m=(mOrder+1)/2;
     
    519526}
    520527
     528void
     529GLInteg::Print(int lp)
     530{
     531 Integrator::Print(lp);
     532 cout<<"GLInteg order="<<mOrder<<endl;
     533 if(lp>0 && mOrder>0) {
     534   for(int i=0;i<mOrder;i++)
     535     cout<<" ("<<mXPos[i]<<","<<mWeights[i]<<")";
     536   cout<<endl;
     537 }
     538}
  • trunk/SophyaLib/NTools/integ.h

    r2283 r2423  
    4343  virtual Integrator& Func(FUNC);
    4444  virtual Integrator& Func(GeneralFunction*, double* par);
     45
     46  virtual inline  int GetNStep(void) {return mNStep;}
     47  virtual void Print(int lp=0);
    4548
    4649protected:
     
    9295  virtual GLInteg& AddBound(double x);
    9396  virtual GLInteg& SetOrder(int order);
     97  virtual inline  int GetOrder(void) {return mOrder;}
     98  virtual void Print(int lp=0);
    9499protected:
    95100
Note: See TracChangeset for help on using the changeset viewer.