Go to the documentation of this file.00001
00002
00003
00004
00005
00006
00011 #include "fvar.hpp"
00012 #ifdef __TURBOC__
00013 #pragma hdrstop
00014 #endif
00015
00016 #include <stdlib.h>
00017 #include <stdio.h>
00018 #include <math.h>
00019
00020 double dmin(double,double);
00021 double dmax(double, double);
00022
00023 #define USE_BARD_PEN
00024
00029 dvariable boundp(dvariable xx, double fmin, double fmax,
00030 const prevariable& _fpen, const double& s)
00031 {
00032 prevariable& fpen=(prevariable&) _fpen;
00033
00034 dvariable t,y,x;
00035 x=xx/s;
00036 const double l4=log(4.0);
00037 dvariable ss=(sin(x*1.57079632679489661)+1.)/2.;
00038 double diff=fmax-fmin;
00039 t=fmin + diff*ss;
00040 #ifdef USE_BARD_PEN
00041
00042 double pen=.00001/diff;
00043 fpen=fpen-pen*(log(ss+1.e-40)+log((1.0-ss)+1.e-40)+l4);
00044 #else
00045
00046 if (x < -.9999)
00047 {
00048 fpen+=(x+0.9999)*(x+0.9999);
00049 }
00050
00051 if (x > 0.9999)
00052 {
00053 fpen+=(x-0.9999)*(x-0.9999);
00054 }
00055
00056 if (x < -1.)
00057 {
00058 fpen+=1000.*(x+1.)*(x+1.);
00059 }
00060
00061 if (x > 1.)
00062 {
00063 fpen+=1000.*(x-1.)*(x-1.);
00064 }
00065 #endif
00066 return(t);
00067 }
00068
00073 double boundp( double xx, double fmin, double fmax, const double& _fpen,
00074 const double& s)
00075 {
00076 double& fpen = (double&)_fpen;
00077 double t,x;
00078 x=xx/s;
00079 const double l4=log(4.0);
00080 double ss=(sin(x*1.57079632679489661)+1.)/2.;
00081 double diff=fmax-fmin;
00082 t=fmin + diff*ss;
00083 #ifdef USE_BARD_PEN
00084
00085 double pen=.00001/diff;
00086 fpen-=pen*(log(ss+1.e-40)+log((1.0-ss)+1.e-40)+l4);
00087 #else
00088 if (x < -.9999)
00089 {
00090 fpen+=(x+0.9999)*(x+0.9999);
00091 }
00092
00093 if (x > 0.9999)
00094 {
00095 fpen+=(x-0.9999)*(x-0.9999);
00096 }
00097
00098 if (x < -1)
00099 {
00100 fpen+=1000*(x+1)*(x+1);
00101 }
00102
00103 if (x > 1)
00104 {
00105 fpen+=1000*(x-1)*(x-1);
00106 }
00107 #endif
00108 return(t);
00109 }
00110
00115 double boundpin(double x, double fmin, double fmax,const double& s)
00116 {
00117 double tinv;
00118
00119 if (x < fmin)
00120 {
00121 if (ad_printf)
00122 (*ad_printf)("variable out of bounds in boundpin: variable = %lg", x);
00123 if (ad_printf) (*ad_printf)("; min = %lg", fmin);
00124 if (ad_printf) (*ad_printf)("; max = %lg\n", fmax);
00125
00126 x=dmin(fmin+.001,fmin+.01*(fmax-fmin));
00127 }
00128
00129 if (x > fmax)
00130 {
00131 if (ad_printf)
00132 (*ad_printf)("variable out of bounds in boundpin: variable = %lg", x);
00133 if (ad_printf) (*ad_printf)("; min = %lg", fmin);
00134 if (ad_printf) (*ad_printf)("; max = %lg\n", fmax);
00135
00136 x=dmax(fmax-.001,fmax-.01*(fmax-fmin));
00137 }
00138
00139 tinv=::asin(2.*(x-fmin)/(fmax-fmin)-1.)/1.57079632679489661;
00140 return(s*tinv);
00141 }
00142
00143
00144
00145
00146
00147
00148
00149
00150
00151
00152
00153
00154
00155
00156
00157
00158
00159
00160
00161
00162
00163
00164
00165
00166
00167
00168
00169
00170
00171
00172
00173