68 const int& nrk_phi)
const {
80 cout <<
"Not yet prepared!!!" << endl ;
87 double xi_t0 = xi_i(0) ;
88 double xi_p0 = xi_i(1) ;
89 double xi_l0 = xi_i(2) ;
92 double dp = 2. * M_PI / double(np) / double(nrk_phi) ;
108 double xi_t1, xi_t2, xi_t3, xi_t4, xi_tf ;
109 double xi_p1, xi_p2, xi_p3, xi_p4, xi_pf ;
110 double xi_l1, xi_l2, xi_l3, xi_l4, xi_lf ;
111 double f1, f2, f3, f4 ;
112 double g1, g2, g3, g4 ;
113 double h1, h2, h3, h4 ;
119 for (
int i=0; i<nrk_phi; i++) {
122 f1 = - xi_l0 * rah * confo2.
val_point(rah, M_PI/2., phi0)
123 + 2. * xi_p0 * dlnconfo.
val_point(rah, M_PI/2., phi0) ;
124 g1 = -2. * xi_t0 * dlnconfo.
val_point(rah, M_PI/2., phi0) ;
125 h1 = (1. - 2.*laplnconfo.
val_point(rah, M_PI/2., phi0)) * xi_t0
126 / rah / confo2.
val_point(rah, M_PI/2., phi0) ;
133 f2 = - (xi_l0+0.5*xi_l1) * rah
134 * confo2.
val_point(rah, M_PI/2., phi0+0.5*dp)
135 + 2. * (xi_p0+0.5*xi_p1)
136 * dlnconfo.
val_point(rah, M_PI/2., phi0+0.5*dp) ;
137 g2 = -2. * (xi_t0+0.5*xi_t1)
138 * dlnconfo.
val_point(rah, M_PI/2., phi0+0.5*dp) ;
139 h2 = (1. - 2.*laplnconfo.
val_point(rah, M_PI/2., phi0+0.5*dp))
140 * (xi_t0+0.5*xi_t1) / rah
141 / confo2.
val_point(rah, M_PI/2., phi0+0.5*dp) ;
148 f3 = - (xi_l0+0.5*xi_l2) * rah
149 * confo2.
val_point(rah, M_PI/2., phi0+0.5*dp)
150 + 2. * (xi_p0+0.5*xi_p2)
151 * dlnconfo.
val_point(rah, M_PI/2., phi0+0.5*dp) ;
152 g3 = -2. * (xi_t0+0.5*xi_t2)
153 * dlnconfo.
val_point(rah, M_PI/2., phi0+0.5*dp) ;
154 h3 = (1. - 2.*laplnconfo.
val_point(rah, M_PI/2., phi0+0.5*dp))
155 * (xi_t0+0.5*xi_t2) / rah
156 / confo2.
val_point(rah, M_PI/2., phi0+0.5*dp) ;
163 f4 = - (xi_l0+xi_l3) * rah * confo2.
val_point(rah, M_PI/2., phi0+dp)
164 + 2. * (xi_p0+xi_p3) * dlnconfo.
val_point(rah, M_PI/2., phi0+dp) ;
165 g4 = -2. * (xi_t0+xi_t3) * dlnconfo.
val_point(rah, M_PI/2., phi0+dp) ;
166 h4 = (1. - 2.*laplnconfo.
val_point(rah, M_PI/2., phi0+dp))
167 * (xi_t0+xi_t3) / rah / confo2.
val_point(rah, M_PI/2., phi0+dp) ;
175 xi_tf = xi_t0 + (xi_t1 + 2.*xi_t2 + 2.*xi_t3 + xi_t4) / 6. ;
176 xi_pf = xi_p0 + (xi_p1 + 2.*xi_p2 + 2.*xi_p3 + xi_p4) / 6. ;
177 xi_lf = xi_l0 + (xi_l1 + 2.*xi_l2 + 2.*xi_l3 + xi_l4) / 6. ;
188 xi_f.
set(0) = xi_tf ;
189 xi_f.
set(1) = xi_pf ;
190 xi_f.
set(2) = xi_lf ;