00001
00002
00003
00004
00005
00006
00007
00008
00009 #ifndef AMROC_EVALUATOR_H
00010 #define AMROC_EVALUATOR_H
00011
00012 #include <iostream>
00013 #include <cstdio>
00014 #include <cmath>
00015 #include <cfloat>
00016
00024 #include "EvaluatorBase.h"
00025
00033 template <class VectorType>
00034 class Evaluator : public EvaluatorBase {
00035 typedef typename VectorType::InternalDataType DataType;
00036 typedef EvaluatorBase base;
00037
00038 public:
00039 typedef GridData<VectorType,3> vector_block_type;
00040 typedef GridData<DataType,3> data_block_type;
00041
00042 Evaluator() : EvaluatorBase() {
00043 std::sprintf(title,"3D-Euler equations");
00044 }
00045
00046 virtual void SetUp(VisualGridBase *gr) {
00047 EvaluatorBase::SetUp(gr);
00048
00049 std::sprintf(&(tkeys[ 0*LEN_TKEYS]), "Distribution ");
00050 std::sprintf(&(tkeys[ 1*LEN_TKEYS]), "Levels ");
00051 std::sprintf(&(tkeys[ 2*LEN_TKEYS]), "Density ");
00052 std::sprintf(&(tkeys[ 3*LEN_TKEYS]), "Velocity u ");
00053 std::sprintf(&(tkeys[ 4*LEN_TKEYS]), "Velocity v ");
00054 std::sprintf(&(tkeys[ 5*LEN_TKEYS]), "Velocity w ");
00055 std::sprintf(&(tkeys[ 6*LEN_TKEYS]), "Total Energy ");
00056 std::sprintf(&(tkeys[ 7*LEN_TKEYS]), "|Velocity| ");
00057 std::sprintf(&(tkeys[ 8*LEN_TKEYS]), "Flow Vectors ");
00058 std::sprintf(&(tkeys[ 9*LEN_TKEYS]), "Schl.-Plot Density ");
00059 std::sprintf(&(tkeys[10*LEN_TKEYS]), "Schl.-Plot Velocity u ");
00060 std::sprintf(&(tkeys[11*LEN_TKEYS]), "Schl.-Plot Velocity v ");
00061 std::sprintf(&(tkeys[12*LEN_TKEYS]), "Schl.-Plot Velocity w ");
00062 std::sprintf(&(tkeys[13*LEN_TKEYS]), "Schl.-Plot Total Energy ");
00063 std::sprintf(&(tkeys[14*LEN_TKEYS]), "Schl.-Plot |Velocity| ");
00064 std::sprintf(&(tkeys[15*LEN_TKEYS]), "Vorticity of w and v ");
00065 std::sprintf(&(tkeys[16*LEN_TKEYS]), "Vorticity of u and w ");
00066 std::sprintf(&(tkeys[17*LEN_TKEYS]), "Vorticity of v and u ");
00067
00068 base::fkeys[0] = base::Grid->FKeys();
00069 base::fkeys[1] = base::Grid->FKeys();
00070 base::fkeys[2] = base::Grid->FKeys();
00071 base::fkeys[3] = base::Grid->FKeys();
00072 base::fkeys[4] = base::Grid->FKeys();
00073 base::fkeys[5] = base::Grid->FKeys();
00074 base::fkeys[6] = base::Grid->FKeys();
00075 base::fkeys[7] = base::Grid->FKeys();
00076 base::fkeys[8] = 2;
00077 base::fkeys[9] = base::Grid->FKeys();
00078 base::fkeys[10] = base::Grid->FKeys();
00079 base::fkeys[11] = base::Grid->FKeys();
00080 base::fkeys[12] = base::Grid->FKeys();
00081 base::fkeys[13] = base::Grid->FKeys();
00082 base::fkeys[14] = base::Grid->FKeys();
00083 base::fkeys[15] = base::Grid->FKeys();
00084 base::fkeys[16] = base::Grid->FKeys();
00085 base::fkeys[17] = base::Grid->FKeys();
00086
00087 base::ikeys[0] = 105; base::ikeys[1] = 115; base::ikeys[2] = 100; base::ikeys[3] = 117;
00088 base::ikeys[4] = 118; base::ikeys[5] = 119; base::ikeys[6] = 101; base::ikeys[7] = 108;
00089 base::ikeys[8] = 102; base::ikeys[9] = 68; base::ikeys[10] = 85; base::ikeys[11] = 86;
00090 base::ikeys[12] = 87; base::ikeys[13] = 69; base::ikeys[14] = 76; base::ikeys[15] = 120;
00091 base::ikeys[16] = 121; base::ikeys[17] = 122;
00092 }
00093
00094 virtual int NKeys() const { return 18; }
00095
00096 virtual void SetScalCells(int key,float v[],int& offset,
00097 vector_block_type& DB, float* dx, bool* CompRead) {
00098 int idx=0;
00099 int keym1 = key-1;
00100 int ac;
00101 BeginFastIndex3(dat, DB.bbox(), DB.data(), VectorType);
00102 switch (key) {
00103 case 3:
00104 case 4:
00105 case 5:
00106 case 6:
00107 case 7:
00108 case 10:
00109 case 11:
00110 case 12:
00111 case 13:
00112 case 14:
00113 ac = ((key>=10 && key<=14) ? 7 : 0);
00114 if (!CompRead[key-(3+ac)])
00115 return;
00116 for_3 (i,j,k,DB.bbox(),DB.bbox().stepsize())
00117 v[idx+offset]= FastIndex3(dat,i,j,k)(key-(3+ac));
00118 mvals[2*keym1] = Min(mvals[2*keym1], v[idx+offset]);
00119 mvals[2*keym1+1] = Max(mvals[2*keym1+1], v[idx+offset]);
00120 idx++;
00121 end_for
00122 break;
00123
00124 case 8:
00125 case 15:
00126 case 16:
00127 case 17:
00128 case 18:
00129 float u1, u2, u3;
00130 for_3 (i,j,k,DB.bbox(),DB.bbox().stepsize())
00131 u1=FastIndex3(dat,i,j,k)(1);
00132 u2=FastIndex3(dat,i,j,k)(2);
00133 u3=FastIndex3(dat,i,j,k)(3);
00134 v[idx+offset]=std::sqrt(u1*u1+u2*u2+u3*u3);
00135 mvals[2*keym1] = Min(mvals[2*keym1], v[idx+offset]);
00136 mvals[2*keym1+1] = Max(mvals[2*keym1+1], v[idx+offset]);
00137 idx++;
00138 end_for
00139 break;
00140
00141 default:
00142 break;
00143 }
00144 EndFastIndex3(dat);
00145 offset+=idx;
00146
00147 if ((key>=10 && key<=14) || (key>=15 && key<=18))
00148 std::cerr << &(tkeys[keym1*LEN_TKEYS])
00149 << " is NOT supported for Type=6. Use Type=1 instead!"
00150 << std::endl;
00151 }
00152
00153 virtual void SetScalNodes(int key,float v[],int& offset,
00154 vector_block_type& DB, const BBox& bbox,
00155 float* dx, bool* CompRead) {
00156 int idx=0;
00157 int keym1 = key-1;
00158 BeginFastIndex3(dat, DB.bbox(), DB.data(), VectorType);
00159 const Coords& step = bbox.stepsize();
00160 switch (key) {
00161 case 3:
00162 case 4:
00163 case 5:
00164 case 6:
00165 case 7:
00166 {
00167 if (!CompRead[key-3])
00168 return;
00169 for_3 (i,j,k,bbox,step)
00170 BBox AddOver(3,i-step(0),j-step(1),k-step(2),
00171 i,j,k,step(0),step(1),step(2));
00172 AddOver *= DB.bbox();
00173
00174 float c=0;
00175 v[idx+offset] = 0.0;
00176 for_3 (l, m, n, AddOver, step)
00177 if (FastIndex3(dat,l,m,n)(key-3) < FLT_MAX) {
00178 v[idx+offset] += FastIndex3(dat,l,m,n)(key-3);
00179 c += 1.0;
00180 }
00181 end_for
00182 if (c>0) v[idx+offset] /= c;
00183
00184 mvals[2*keym1] = Min(mvals[2*keym1], v[idx+offset]);
00185 mvals[2*keym1+1] = Max(mvals[2*keym1+1], v[idx+offset]);
00186 idx++;
00187 end_for
00188 break;
00189 }
00190 case 10:
00191 case 11:
00192 case 12:
00193 case 13:
00194 case 14:
00195 {
00196 if (!CompRead[key-10])
00197 return;
00198 for_3 (i,j,k,bbox,step)
00199 BBox AddOver(3,i-step(0),j-step(1),k-step(2),
00200 i,j,k,step(0),step(1),step(2));
00201 AddOver *= DB.bbox();
00202
00203 data_block_type DBHelp(AddOver);
00204 BeginFastIndex3(help, DBHelp.bbox(), DBHelp.data(), DataType);
00205 for_3 (l, m, n, AddOver, step)
00206 FastIndex3(help,l,m,n) = FastIndex3(dat,l,m,n)(key-10);
00207 end_for
00208 EndFastIndex3(help);
00209 v[idx+offset] = SetScalGradNodes(DBHelp,dx);
00210
00211 mvals[2*keym1] = Min(mvals[2*keym1], v[idx+offset]);
00212 mvals[2*keym1+1] = Max(mvals[2*keym1+1], v[idx+offset]);
00213 idx++;
00214 end_for
00215 break;
00216 }
00217 case 8:
00218 {
00219 float u1, u2, u3;
00220 for_3 (i, j, k, bbox, step)
00221 BBox AddOver(3,i-step(0),j-step(1),k-step(2),
00222 i,j,k,step(0),step(1),step(2));
00223 AddOver *= DB.bbox();
00224
00225 float c=0;
00226 v[idx+offset] = 0.0;
00227 for_3 (l, m, n, AddOver, step)
00228 u1 = FastIndex3(dat,l,m,n)(1);
00229 u2 = FastIndex3(dat,l,m,n)(2);
00230 u3 = FastIndex3(dat,l,m,n)(3);
00231 if (u1<FLT_MAX && u2<FLT_MAX && u3<FLT_MAX) {
00232 v[idx+offset] += std::sqrt(u1*u1+u2*u2+u3*u3);
00233 c += 1.0;
00234 }
00235 end_for
00236 if (c>0) v[idx+offset] /= c;
00237
00238 mvals[2*keym1] = Min(mvals[2*keym1], v[idx+offset]);
00239 mvals[2*keym1+1] = Max(mvals[2*keym1+1], v[idx+offset]);
00240 idx++;
00241 end_for
00242 break;
00243 }
00244 case 15:
00245 {
00246 float u1, u2, u3;
00247 for_3 (i, j, k, bbox, step)
00248 BBox AddOver(3,i-step(0),j-step(1),k-step(2),
00249 i,j,k,step(0),step(1),step(2));
00250 AddOver *= DB.bbox();
00251
00252 data_block_type DBHelp(AddOver);
00253 DBHelp = FLT_MAX;
00254 BeginFastIndex3(help, DBHelp.bbox(), DBHelp.data(), DataType);
00255 for_3 (l, m, n, AddOver, step)
00256 u1 = FastIndex3(dat,l,m,n)(1);
00257 u2 = FastIndex3(dat,l,m,n)(2);
00258 u3 = FastIndex3(dat,l,m,n)(3);
00259 if (u1<FLT_MAX && u2<FLT_MAX && u3<FLT_MAX)
00260 FastIndex3(help,l,m,n) = std::sqrt(u1*u1+u2*u2+u3*u3);
00261 end_for
00262 EndFastIndex3(help);
00263 v[idx+offset] = SetScalGradNodes(DBHelp,dx);
00264
00265 mvals[2*keym1] = Min(mvals[2*keym1], v[idx+offset]);
00266 mvals[2*keym1+1] = Max(mvals[2*keym1+1], v[idx+offset]);
00267 idx++;
00268 end_for
00269 break;
00270 }
00271 case 16:
00272 {
00273 if (!CompRead[2] || !CompRead[3])
00274 return;
00275 for_3 (i,j,k,bbox,step)
00276 BBox AddOver(3,i-step(0),j-step(1),k-step(2),
00277 i,j,k,step(0),step(1),step(2));
00278 AddOver *= DB.bbox();
00279
00280 data_block_type u3(AddOver), u2(AddOver);
00281 BeginFastIndex3(u2, u2.bbox(), u2.data(), DataType);
00282 BeginFastIndex3(u3, u3.bbox(), u3.data(), DataType);
00283 for_3 (l, m, n, AddOver, step)
00284 FastIndex3(u2,l,m,n) = FastIndex3(dat,l,m,n)(2);
00285 FastIndex3(u3,l,m,n) = FastIndex3(dat,l,m,n)(3);
00286 end_for
00287 EndFastIndex3(u2);
00288 EndFastIndex3(u3);
00289 v[idx+offset] = ScalDerivNodes(u3,dx,1)-ScalDerivNodes(u2,dx,2);
00290
00291 mvals[2*keym1] = Min(mvals[2*keym1], v[idx+offset]);
00292 mvals[2*keym1+1] = Max(mvals[2*keym1+1], v[idx+offset]);
00293 idx++;
00294 end_for
00295 break;
00296 }
00297 case 17:
00298 {
00299 if (!CompRead[1] || !CompRead[3])
00300 return;
00301 for_3 (i,j,k,bbox,step)
00302 BBox AddOver(3,i-step(0),j-step(1),k-step(2),
00303 i,j,k,step(0),step(1),step(2));
00304 AddOver *= DB.bbox();
00305
00306 data_block_type u1(AddOver), u3(AddOver);
00307 BeginFastIndex3(u1, u1.bbox(), u1.data(), DataType);
00308 BeginFastIndex3(u3, u3.bbox(), u3.data(), DataType);
00309 for_3 (l, m, n, AddOver, step)
00310 FastIndex3(u1,l,m,n) = FastIndex3(dat,l,m,n)(1);
00311 FastIndex3(u3,l,m,n) = FastIndex3(dat,l,m,n)(3);
00312 end_for
00313 EndFastIndex3(u1);
00314 EndFastIndex3(u3);
00315 v[idx+offset] = ScalDerivNodes(u1,dx,2)-ScalDerivNodes(u3,dx,0);
00316
00317 mvals[2*keym1] = Min(mvals[2*keym1], v[idx+offset]);
00318 mvals[2*keym1+1] = Max(mvals[2*keym1+1], v[idx+offset]);
00319 idx++;
00320 end_for
00321 break;
00322 }
00323 case 18:
00324 {
00325 if (!CompRead[1] || !CompRead[2])
00326 return;
00327 for_3 (i,j,k,bbox,step)
00328 BBox AddOver(3,i-step(0),j-step(1),k-step(2),
00329 i,j,k,step(0),step(1),step(2));
00330 AddOver *= DB.bbox();
00331
00332 data_block_type u1(AddOver), u2(AddOver);
00333 BeginFastIndex3(u1, u1.bbox(), u1.data(), DataType);
00334 BeginFastIndex3(u2, u2.bbox(), u2.data(), DataType);
00335 for_3 (l, m, n, AddOver, step)
00336 FastIndex3(u1,l,m,n) = FastIndex3(dat,l,m,n)(1);
00337 FastIndex3(u2,l,m,n) = FastIndex3(dat,l,m,n)(2);
00338 end_for
00339 EndFastIndex3(u1);
00340 EndFastIndex3(u2);
00341 v[idx+offset] = ScalDerivNodes(u2,dx,0)-ScalDerivNodes(u1,dx,1);
00342
00343 mvals[2*keym1] = Min(mvals[2*keym1], v[idx+offset]);
00344 mvals[2*keym1+1] = Max(mvals[2*keym1+1], v[idx+offset]);
00345 idx++;
00346 end_for
00347 break;
00348 }
00349
00350 default:
00351 break;
00352 }
00353 EndFastIndex3(dat);
00354 offset+=idx;
00355 }
00356
00357 protected:
00358 DataType ScalDerivNodes(data_block_type& DB, float* dx, const int d) {
00359 float grad = 0.0;
00360 float c=0.0;
00361 BeginFastIndex3(dat, DB.bbox(), DB.data(), DataType);
00362 const BBox& AddOver = DB.bbox();
00363 float wdx[3];
00364 float v1, v2;
00365 for (int i=0; i<3; i++)
00366 wdx[i] = dx[i]*AddOver.stepsize(i);
00367
00368 switch(d) {
00369 case 0:
00370 if (AddOver.extents(0) == 2) {
00371 for (int m=AddOver.lower(1); m<=AddOver.upper(1);
00372 m+=AddOver.stepsize(1)) {
00373 for (int n=AddOver.lower(2); n<=AddOver.upper(2);
00374 n+=AddOver.stepsize(2)) {
00375 v1 = FastIndex3(dat,AddOver.lower(0),m,n);
00376 v2 = FastIndex3(dat,AddOver.upper(0),m,n);
00377 if (v1<FLT_MAX && v2<FLT_MAX) {
00378 grad += (v2-v1)/wdx[0];
00379 c += 1.0;
00380 }
00381 }
00382 }
00383 if (c>0.0) grad /= c;
00384 }
00385 break;
00386 case 1:
00387 if (AddOver.extents(1) == 2) {
00388 for (int l=AddOver.lower(0); l<=AddOver.upper(0);
00389 l+=AddOver.stepsize(0)) {
00390 for (int n=AddOver.lower(2); n<=AddOver.upper(2);
00391 n+=AddOver.stepsize(2)) {
00392 v1 = FastIndex3(dat,l,AddOver.lower(1),n);
00393 v2 = FastIndex3(dat,l,AddOver.upper(1),n);
00394 if (v1<FLT_MAX && v2<FLT_MAX) {
00395 grad += (v2-v1)/wdx[1];
00396 c += 1.0;
00397 }
00398 }
00399 }
00400 if (c>0.0) grad /= c;
00401 }
00402 break;
00403 case 2:
00404 if (AddOver.extents(2) == 2) {
00405 for (int l=AddOver.lower(0); l<=AddOver.upper(0);
00406 l+=AddOver.stepsize(0)) {
00407 for (int m=AddOver.lower(1); m<=AddOver.upper(1);
00408 m+=AddOver.stepsize(1)) {
00409 v1 = FastIndex3(dat,l,m,AddOver.lower(2));
00410 v2 = FastIndex3(dat,l,m,AddOver.upper(2));
00411 if (v1<FLT_MAX && v2<FLT_MAX) {
00412 grad += (v2-v1)/wdx[2];
00413 c += 1.0;
00414 }
00415 }
00416 }
00417 if (c>0.0) grad /= c;
00418 }
00419 break;
00420 default:
00421 break;
00422 }
00423 EndFastIndex3(dat);
00424
00425 return grad;
00426 }
00427
00428 DataType SetScalGradNodes(data_block_type& DB, float* dx) {
00429 float xgrad, ygrad, zgrad;
00430 xgrad = ScalDerivNodes(DB,dx,0);
00431 ygrad = ScalDerivNodes(DB,dx,1);
00432 zgrad = ScalDerivNodes(DB,dx,2);
00433 return std::sqrt(xgrad*xgrad+ygrad*ygrad+zgrad*zgrad);
00434 }
00435 };
00436
00437
00438 #endif
