Difference between revisions of "LogZp1oZm1"

From aemwiki
Jump to: navigation, search
m
(revised summary of speedup, changed y=0 to y>=0 for consistency with c++ complex library)
Line 3: Line 3:
 
; Purpose: Calculates <math>\ln\left(\frac{Z-1}{Z+1}\right)</math>
 
; Purpose: Calculates <math>\ln\left(\frac{Z-1}{Z+1}\right)</math>
  
; Author: James R. Craig
+
; Author: James R. Craig, Ph. Le Grand
  
 
  #include <math.h>
 
  #include <math.h>
Line 12: Line 12:
 
  logZm1oZp1
 
  logZm1oZp1
 
  optimized version of log ((Z-1.0)/(Z+1.0)):  
 
  optimized version of log ((Z-1.0)/(Z+1.0)):  
  1.1% of the cost of a direct call- 2 orders of magnitude speedup!
+
  ~13% of the cost of a direct call- 1 order of magnitude speedup!
 
  Author: James R. Craig
 
  Author: James R. Craig
 
  Modified 2007 12 06 by Ph. Le Grand:
 
  Modified 2007 12 06 by Ph. Le Grand:
Line 28: Line 28:
 
  if ((r2)<1.0){
 
  if ((r2)<1.0){
 
  add=-PI;
 
  add=-PI;
  if (y>0.0){add=PI;}
+
  if (y>=0.0){add=PI;}
 
  }
 
  }
 
  return cmplex(0.5*log(1.0-((x*(den-2)+den-r2-1.0)/(0.25*den*den))),
 
  return cmplex(0.5*log(1.0-((x*(den-2)+den-r2-1.0)/(0.25*den*den))),
 
        add+atan( (y+y) /(den-x-x-2.0) ));   
 
        add+atan( (y+y) /(den-x-x-2.0) ));   
 
  }
 
  }

Revision as of 08:51, 7 December 2007

Language
c++
Purpose
Calculates <math>\ln\left(\frac{Z-1}{Z+1}\right)</math>
Author
James R. Craig, Ph. Le Grand
#include <math.h>
#include <complex>
typedef std::complex<double> cmplex;
const double PI=3.141592653589793238462643; // pi 
/*******************************************************************************
logZm1oZp1
optimized version of log ((Z-1.0)/(Z+1.0)): 
~13% of the cost of a direct call- 1 order of magnitude speedup!
Author: James R. Craig
Modified 2007 12 06 by Ph. Le Grand:
    Added local variables x,y,r2 to reduce number of calls to 
    member functions of Z and reduce number of operations 
    for minor speed gains.
-------------------------------------------------------------------------------*/
inline cmplex logZm1oZp1(const cmplex &Z)
{ 	
	double x(Z.real());
	double y(Z.imag());	
	double r2=x*x+y*y;								 
	double den(r2+x+x+1.0);
	double add=0.0; 
	if ((r2)<1.0){
		add=-PI;
		if (y>=0.0){add=PI;}
	}
	return cmplex(0.5*log(1.0-((x*(den-2)+den-r2-1.0)/(0.25*den*den))),
		      add+atan( (y+y) /(den-x-x-2.0) ));   
}