CirclePattern.h

Go to the documentation of this file.
00001 
00005 #ifndef CIRCLE_PATTERN_CM
00006 #define CIRCLE_PATTERN_CM
00007 #define NOMINMAX
00008 #include <iostream>
00009 #define _USE_MATH_DEFINES
00010 #include <math.h>
00011 
00012 class Mesh;
00013 
00043 struct CirclePattern {
00044 
00045    static double evaluateEnergy (Mesh * mesh, double * xx);
00046 
00047    static void computeGradient (Mesh * mesh, double * xx, double * grad);
00048 
00049    static void buildSparsityForHessian (Mesh * mesh, int * hessI, int * hessJ);
00050 
00051    static void computeHessian (Mesh * mesh, double * xx, double * hessval);
00052 
00053 
00054    static double ImLi2Sum (const double& theta, const double& x, 
00055       const double& tanHalfThetaStar, const double& clThetaStar);
00056 
00058    static double Cl2 (const double & _x);
00059 
00062    static double fTheta2(const double& cosTheta, const double& sinTheta, const double& x);
00063 
00065    static double fThetaPrime2(const double & cosTheta, const double& sinTheta, const double & x) {
00066       double res = sinTheta/(cosh(x) - cosTheta);
00067       return res;
00068    }
00069 
00071    static double IEEEremainder(const double & x, const double & div);
00072 };
00073 
00074 inline double CirclePattern::fTheta2(const double& cosTheta, const double& sinTheta, const double& x){
00075    double ex = exp(x);
00076    double s1 = 1.0 - ex*cosTheta;
00077    double s2 = ex*sinTheta;
00078 
00079    double res = M_PI - 2*atan2(s1, s2);
00080 
00081    return res;
00082 }
00083 
00084 inline double CirclePattern::IEEEremainder(const double & x, const double & div) {
00085 
00086    double  res;
00087    if (x >= 0) {
00088       if (x < div) res = x;
00089       else {
00090           double part = floor(x / div);
00091           res = x - part*div;
00092       }
00093    }
00094    else {
00095       double part = floor(x / div);
00096       res = x - part*div;
00097    }
00098    if (res < 0 || res > div) {
00099       std::cout << "x " << x << " div " << div << std::endl;   
00100       std::cerr << "Wrong IEEEreminder!!!" << std::endl;
00101    }
00102    return res;
00103 }
00104 
00105 inline double CirclePattern::Cl2 (const double & _x) {
00106    double TWO_PI = 2 * M_PI;
00107    double MINUS_LOG_2 = log(0.5);
00108    double BORDERLINE = 2.0944;
00109 
00110    if (_x == 0.0) {
00111       return 0.0;
00112    }
00113 
00114    double x = IEEEremainder(_x, TWO_PI);
00115 
00116    if (x == 0.0) {
00117       return 0.0;
00118    }
00119 
00120    if (fabs(x) <= BORDERLINE) {
00121       double xx = x * x;
00122       return ((((((((((((
00123           2.3257441143020875e-22 * xx
00124           + 1.0887357368300848e-20) * xx 
00125           + 5.178258806090624e-19) * xx
00126           + 2.5105444608999545e-17) * xx
00127           + 1.2462059912950672e-15) * xx
00128           + 6.372636443183181e-14) * xx
00129           + 3.387301370953521e-12) * xx
00130           + 1.8978869988971e-10) * xx
00131           + 1.1482216343327455e-8) * xx
00132           + 7.873519778281683e-7) * xx
00133           + 0.00006944444444444444) * xx
00134           + 0.013888888888888888) * xx 
00135           - log(fabs(x)) + 1.0) * x;
00136    }
00137    x += ((x > 0.0) ? - M_PI : M_PI);
00138    double xx = x * x;
00139    return ((((((((((((  
00140       3.901950904063069e-15 * xx
00141       + 4.566487567193635e-14) * xx
00142       + 5.429792727596476e-13) * xx
00143       + 6.5812165661369675e-12) * xx
00144       + 8.167010963952222e-11) * xx
00145       + 1.0440290284867003e-9) * xx
00146       + 1.3870999114054669e-8) * xx
00147       + 1.941538399871733e-7) * xx
00148       + 2.927965167548501e-6) * xx
00149       + 0.0000496031746031746) * xx
00150       + 0.0010416666666666667) * xx
00151       + 0.041666666666666664) * xx 
00152       + MINUS_LOG_2) * x;
00153 }
00154 
00155 #endif
00156 

Generated on Sat Jun 3 13:33:41 2006 for CirclePatterns by  doxygen 1.4.5