Changeset 921 for trunk/source/global/management/src/G4PhysicsVector.cc

Timestamp:

Feb 16, 2009, 10:14:30 AM (15 years ago)

Author:

garnier

Message:

en test de gl2ps. Problemes de libraries

File:

• trunk/source/global/management/src/G4PhysicsVector.cc (modified) (7 diffs)

trunk/source/global/management/src/G4PhysicsVector.cc

@@ -25 +25 @@
 //
 //
-// $Id: G4PhysicsVector.cc,v 1.22 2008/09/05 18:04:45 vnivanch Exp $
-// GEANT4 tag $Name: HEAD $
+// $Id: G4PhysicsVector.cc,v 1.27 2008/10/16 12:14:36 gcosmo Exp $
+// GEANT4 tag $Name: geant4-09-02-cand-01 $
 //
 // 
@@ -52 +52 @@
 G4PhysicsVector::G4PhysicsVector(G4bool spline)
  : type(T_G4PhysicsVector),
-   edgeMin(DBL_MAX), edgeMax(0.), numberOfBin(0),
+   edgeMin(0.), edgeMax(0.), numberOfBin(0),
    lastEnergy(0.), lastValue(0.), lastBin(0), 
    secDerivative(0), useSpline(spline)
-{
-  binVector.push_back(edgeMin);
-  dataVector.push_back(0.0);
-}
+{}
 
 // --------------------------------------------------------------
@@ -64 +61 @@
 G4PhysicsVector::~G4PhysicsVector() 
 {
-  delete [] secDerivative;
+  DeleteData();
 }
 
@@ -71 +68 @@
 G4PhysicsVector::G4PhysicsVector(const G4PhysicsVector& right)
 {
-  *this=right;
+  CopyData(right);
 }
 
@@ -81 +78 @@
   if (type != right.type)  { return *this; }
 
-  type = right.type;
-  edgeMin = right.edgeMin;
-  edgeMax = right.edgeMax;
-  numberOfBin = right.numberOfBin;
-  lastEnergy = right.lastEnergy;
-  lastValue = right.lastValue;
-  lastBin = right.lastBin;
-  dataVector = right.dataVector;
-  binVector = right.binVector;
-  secDerivative = right.secDerivative;
-  useSpline = right.useSpline;
-  comment = right.comment;
+  DeleteData();
+  CopyData(right);
+
   return *this;
 }
@@ -108 +96 @@
 {
   return (this != &right);
+}
+
+// --------------------------------------------------------------
+
+void G4PhysicsVector::DeleteData()
+{
+  delete [] secDerivative;
+  secDerivative = 0;
+}
+
+// --------------------------------------------------------------
+
+void G4PhysicsVector::CopyData(const G4PhysicsVector& vec)
+{
+  type = vec.type;
+  edgeMin = vec.edgeMin;
+  edgeMax = vec.edgeMax;
+  numberOfBin = vec.numberOfBin;
+  lastEnergy = vec.lastEnergy;
+  lastValue = vec.lastValue;
+  lastBin = vec.lastBin;
+  dataVector = vec.dataVector;
+  binVector = vec.binVector;
+  useSpline = vec.useSpline;
+  comment = vec.comment;
+  if (vec.secDerivative)
+  {
+    secDerivative = new G4double [numberOfBin];
+    for (size_t i=0; i<numberOfBin; i++)
+    {
+       secDerivative[i] = vec.secDerivative[i];
+    }
+  }
+  else
+  {
+    secDerivative = 0;
+  }
 }
 
@@ -224 +249 @@
   secDerivative = new G4double [numberOfBin];
 
-  secDerivative[0] = 0.0 ;
+  size_t n = numberOfBin-1;
 
   // cannot compute derivatives for less than 3 points
-  if(3 > numberOfBin) {
-    secDerivative[numberOfBin-1] = 0.0 ;
+  if(3 > numberOfBin)
+  {
+    secDerivative[0] = 0.0;
+    secDerivative[n] = 0.0;
     return;
   }
 
-  G4double* u = new G4double [numberOfBin];
-  u[0] = 0.0 ;
-   
-  // Decomposition loop for tridiagonal algorithm. secDerivative[i]
-  // and u[i] are used for temporary storage of the decomposed factors.
-   
-  for(size_t i=1; i<numberOfBin-1; i++)
-  {
-    G4double sig = (binVector[i]-binVector[i-1])
-                 / (binVector[i+1]-binVector[i-1]) ;
-    G4double p = sig*secDerivative[i-1] + 2.0 ;
-    secDerivative[i] = (sig - 1.0)/p ;
-    u[i] = (dataVector[i+1]-dataVector[i])/(binVector[i+1]-binVector[i])
-         - (dataVector[i]-dataVector[i-1])/(binVector[i]-binVector[i-1]) ;
-    u[i] =(6.0*u[i]/(binVector[i+1]-binVector[i-1]) - sig*u[i-1])/p ;
-  }
-
-  G4double qn = 0.0 ;
-  G4double un = 0.0 ;
-
-  secDerivative[numberOfBin-1] = (un - qn*u[numberOfBin-2])
-                               / (qn*secDerivative[numberOfBin-2] + 1.0) ;
-   
-  // The back-substitution loop for the triagonal algorithm of solving
-  // a linear system of equations.
-   
-  for(G4int k=numberOfBin-2; k>=0; k--)
-  {
-    secDerivative[k] = secDerivative[k]*secDerivative[k+1] + u[k];
-  }
-  delete [] u;
+  for(size_t i=1; i<n; i++)
+  {
+    secDerivative[i] = 
+      3.0*((dataVector[i+1]-dataVector[i])/(binVector[i+1]-binVector[i]) -
+           (dataVector[i]-dataVector[i-1])/(binVector[i]-binVector[i-1]))
+      /(binVector[i+1]-binVector[i-1]);
+  }
+  secDerivative[n] = secDerivative[n-1];
+  secDerivative[0] = secDerivative[1];
 }