C++ std::imag方法代碼示例

// integral of product of exponential terms arising in Green's function
// computations (i.e. integral of term in zsProd() from 0 to ds)
cdouble zsInt(int s1, int s2, cdouble kzs, double ds) {
	if(s1==1 && s2==1)
		return (exp(2. * II * real(kzs) * ds) - 1.) / (2. * II * real(kzs));
	else if(s1==1 && s2==-1)
		return (1. - exp(-2. * imag(kzs) * ds)) / (2. * imag(kzs));
	else if(s1==-1 && s2==1)
		return (exp(2. * imag(kzs) * ds) - 1.) / (2. * imag(kzs));
	else if(s1==-1 && s2==-1)
		return (1. - exp(-2. * II * real(kzs) * ds)) / (2. * II * real(kzs));
	return -1; // shouldn't get here
}

double flux(const mlgeo &g, const SMatrix &S, int l, double zl, 
			const int *s, int Ns, double nHat) { 
	pwaves pTE, pTM;
	double f = 0;
	for (int sind = 0, si = s[sind]; sind < Ns; ++sind, ++si) {
		if (imag(g.eps(si)) == 0)
			continue;
		pWavesL(S, l, si, TE, &pTE);
		pWavesL(S, l, si, TM, &pTM);
		f -= nHat * S.k0 * S.k0 / (M_PI * M_PI) * imag(g.eps(si)) 
			* imag( gfFluxTE(S, pTE, l, si, zl, g.d(si), nHat)
			+ gfFluxTM(S, pTM, l, si, zl, g.d(si), nHat) );
	}
	return f;
}

// same as above but with S precomputed (e.g. if multiple zl's desired)
double dos(const SMatrix &S, int l, double zl) {
	pwaves pTE, pTM;	
	pWavesL(S, l, l, TE, &pTE); // source layer = emitter layer
	pWavesL(S, l, l, TM, &pTM);

	cdouble kzl = S.kz[l]; // kzs = kzl
	cdouble kl = S.k[l];
	cdouble kp = S.kp;
	cdouble xl = II * kzl * zl; // xs = xl
	
	cdouble Ae = pTE.Al;
	cdouble Be = pTE.Bl * exp(-2.*xl);
	cdouble Ce = pTE.Cl * exp(2.*xl);
	cdouble De = pTE.Dl;
	cdouble Am = pTM.Al;
	cdouble Bm = pTM.Bl * exp(-2.*xl);
	cdouble Cm = pTM.Cl * exp(2.*xl);
	cdouble Dm = pTM.Dl;
		
	// extra factors of 1 are extra source term when l==s
	cdouble Epp = II * kzl * kp / (kl * kl) 
				* (Am - Bm - Cm + Dm + 1.); 
	cdouble Ett = II * kp / kzl
				* (Ae + Be + Ce + De + 1.);
	cdouble Ezz = II * kp * kp * kp / (kzl * kl * kl) 
				* (Am + Bm + Cm + Dm + 1.);
	return S.k0 * S.k0 * imag(Epp + Ett + Ezz) / (2. * M_PI * M_PI);   
}

  /** Kernel L2T operation
   * r += Op(L, t) where L is the local expansion and r is the result
   *
   * @param[in] L The local expansion
   * @param[in] center The center of the box with the local expansion
   * @param[in] target The target of this L2T operation
   * @param[in] result The result to accumulate into
   * @pre L includes the influence of all sources outside its box
   */
  void L2T(const local_type& L, const point_type& center,
           const target_type& target, result_type& result) const {
    using std::real;
    using std::imag;

    real_type rho, theta, phi;
    SphOp::cart2sph(rho, theta, phi, target - center);
    complex_type Z[P*(P+1)/2], dZ[P*(P+1)/2];
    SphOp::evalZ(rho, theta, phi, P, Z, dZ);

    point_type sph = point_type();
    int nm = 0;
    for (int n = 0; n != P; ++n) {
      const real_type LZ = real(L[nm])*real(Z[nm]) - imag(L[nm])*imag(Z[nm]);
      result[0] += LZ;
      sph[0]    += LZ / rho * n;
      sph[1]    += real(L[nm])*real(dZ[nm]) - imag(L[nm])*imag(dZ[nm]);

      ++nm;
      for (int m = 1; m <= n; ++m, ++nm) {
        const complex_type LZ = L[nm] * Z[nm];
        result[0] += 2 * real(LZ);
        sph[0]    += 2 * real(LZ) / rho * n;
        sph[1]    += 2 * (real(L[nm])*real(dZ[nm])-imag(L[nm])*imag(dZ[nm]));
        sph[2]    += 2 *-imag(LZ) * m;
      }
    }

    const point_type cart = SphOp::sph2cart(rho, theta, phi, sph);
    result[1] += cart[0];
    result[2] += cart[1];
    result[3] += cart[2];
  }