?? genericpradosschemesforfastmarching_3d.h
字號:
} // Initialisation of Root. T Root[3][3][3]; for(s1=0; s1<3; s1++) for(s2=0; s2<3; s2++) for(s3=0; s3<3; s3++) { Root[s1][s2][s3] = this->big; } // Initialisation of the optimal dynamics. // In dimension 3, the dynamics (function f associated to an HJB eq.) is a 3D-vector. // optDymamics1, optDymamics2, optDymamics3 represente the its componantes. // As for "Root" and "DoExistSol", we define and compute an optimal // Dymamics on each simplex and each "boundaries" of simplexes. T optDymamics1[3][3][3]; T optDymamics2[3][3][3]; T optDymamics3[3][3][3]; for(s1=0; s1<3; s1++) for(s2=0; s2<3; s2++) for(s3=0; s3<3; s3++) { optDymamics1[s1][s2][s3] = 0; optDymamics2[s1][s2][s3] = 0; optDymamics3[s1][s2][s3] = 0; } // =========================================== // Computation of Roots such that opt control // is in interior of the various Ds1s2s3: for(s1=-1; s1<=1; s1+=2) for(s2=-1; s2<=1; s2+=2) for(s3=-1; s3<=1; s3+=2) { DoExistSol[s1+1][s2+1][s3+1] = eqSolverOnPart_with_s1s2s3_nonNull( this->_GetValue(x+s1,y,z), this->_GetValue(x,y+s2,z), this->_GetValue(x,y,z+s3), s1,s2,s3, x,y,z, Root[s1+1][s2+1][s3+1], optDymamics1[s1+1][s2+1][s3+1],optDymamics2[s1+1][s2+1][s3+1], optDymamics3[s1+1][s2+1][s3+1] ); } // =========================================== // Computation of Roots such that opt control // is in interior of the intersection of Ds1s2s3 // and set of a such that f_i(x,a) = 0 s1=0; for(s2=-1; s2<=1; s2+=2) for(s3=-1; s3<=1; s3+=2) { if ( ! (DoExistSol[-1+1][s2+1][s3+1] || DoExistSol[1+1][s2+1][s3+1]) ) { DoExistSol[s1+1][s2+1][s3+1] = eqSolverOnPart_withOne_si_Null( this->_GetValue(x+s1,y,z), this->_GetValue(x,y+s2,z), this->_GetValue(x,y,z+s3), s1,s2,s3, x,y,z, 1, // Indice_i_Such_siEqualTo0; Root[s1+1][s2+1][s3+1], optDymamics1[s1+1][s2+1][s3+1],optDymamics2[s1+1][s2+1][s3+1], optDymamics3[s1+1][s2+1][s3+1] ); } else{ DoExistSol[s1+1][s2+1][s3+1] = true; // useful for the next step! } } s2=0; for(s1=-1; s1<=1; s1+=2) for(s3=-1; s3<=1; s3+=2) { if ( ! (DoExistSol[s1+1][-1+1][s3+1] || DoExistSol[s1+1][1+1][s3+1]) ) { DoExistSol[s1+1][s2+1][s3+1] = eqSolverOnPart_withOne_si_Null( this->_GetValue(x+s1,y,z), this->_GetValue(x,y+s2,z), this->_GetValue(x,y,z+s3), s1,s2,s3, x,y,z, 2, // Indice_i_Such_siEqualTo0; Root[s1+1][s2+1][s3+1], optDymamics1[s1+1][s2+1][s3+1],optDymamics2[s1+1][s2+1][s3+1], optDymamics3[s1+1][s2+1][s3+1] ); } else{ DoExistSol[s1+1][s2+1][s3+1] = true; // useful for the next step! } } s3=0; for(s2=-1; s2<=1; s2+=2) for(s1=-1; s1<=1; s1+=2) { if ( ! (DoExistSol[s1+1][s2+1][-1+1] || DoExistSol[s1+1][s2+1][1+1]) ) { DoExistSol[s1+1][s2+1][s3+1] = eqSolverOnPart_withOne_si_Null( this->_GetValue(x+s1,y,z), this->_GetValue(x,y+s2,z), this->_GetValue(x,y,z+s3), s1,s2,s3, x,y,z, 3, // Indice_i_Such_siEqualTo0; Root[s1+1][s2+1][s3+1], optDymamics1[s1+1][s2+1][s3+1],optDymamics2[s1+1][s2+1][s3+1], optDymamics3[s1+1][s2+1][s3+1] ); } else { DoExistSol[s1+1][s2+1][s3+1] = true; // useful for the next step! } } // =========================================== // computation of Roots such that opt control // is in Interior of the intersection of Ds1s2s3 // and set of a such that f_i(x,a) = 0 s1=0; s2=0; for(s3=-1; s3<=1; s3+=2) { // More comments could be required here! .... if ( ! ( DoExistSol[ 0+1][-1+1][s3+1] || DoExistSol[ 0+1][ 1+1][s3+1] || DoExistSol[-1+1][ 0+1][s3+1] || DoExistSol[ 1+1][ 0+1][s3+1] ) ) { DoExistSol[s1+1][s2+1][s3+1] = eqSolverOnPart_withTwo_si_Null( this->_GetValue(x+s1,y,z), this->_GetValue(x,y+s2,z), this->_GetValue(x,y,z+s3), s1,s2,s3, x,y,z, 3, // Indice_i_Such_siDifferentTo0; Root[s1+1][s2+1][s3+1], optDymamics1[s1+1][s2+1][s3+1],optDymamics2[s1+1][s2+1][s3+1], optDymamics3[s1+1][s2+1][s3+1] ); } } s1=0; s3=0; for(s2=-1; s2<=1; s2+=2) { // More comments could be required here! .... if ( ! ( DoExistSol[ 0+1][s2+1][-1+1] || DoExistSol[ 0+1][s2+1][ 1+1] || DoExistSol[-1+1][s2+1][ 0+1] || DoExistSol[ 1+1][s2+1][ 0+1] ) ) { DoExistSol[s1+1][s2+1][s3+1] = eqSolverOnPart_withTwo_si_Null( this->_GetValue(x+s1,y,z), this->_GetValue(x,y+s2,z), this->_GetValue(x,y,z+s3), s1,s2,s3, x,y,z, 2, // Indice_i_Such_siDifferentTo0; Root[s1+1][s2+1][s3+1], optDymamics1[s1+1][s2+1][s3+1],optDymamics2[s1+1][s2+1][s3+1], optDymamics3[s1+1][s2+1][s3+1] ); } } s2=0; s3=0; for(s1=-1; s1<=1; s1+=2) { // More comments could be required here! .... if ( ! ( DoExistSol[s1+1][ 0+1][-1+1] || DoExistSol[s1+1][ 0+1][ 1+1] || DoExistSol[s1+1][-1+1][ 0+1] || DoExistSol[s1+1][ 1+1][ 0+1] ) ) { DoExistSol[s1+1][s2+1][s3+1] = eqSolverOnPart_withTwo_si_Null( this->_GetValue(x+s1,y,z), this->_GetValue(x,y+s2,z), this->_GetValue(x,y,z+s3), s1,s2,s3, x,y,z, 1, // Indice_i_Such_siDifferentTo0; Root[s1+1][s2+1][s3+1], optDymamics1[s1+1][s2+1][s3+1],optDymamics2[s1+1][s2+1][s3+1], optDymamics3[s1+1][s2+1][s3+1] ); } } // =========================================================================== // The solution is then the smallest Root: T minRoot = this->big; // the smallest root for(s1=0; s1<3; s1++) for(s2=0; s2<3; s2++) for(s3=0; s3<3; s3++) if ((minRoot > Root[s1][s2][s3])) { minRoot = Root[s1][s2][s3]; // the current min is then Root[s1][s2][s3] gs1=s1; gs2=s2; gs3=s3; // remember the simplex associated to the current min }; // Let us record the optimal Dymamics associated to the adequate simplex. // The packed optimal dynamics is the one associated to the used simplex. // Note: Preservation of the optimal Dymamics is not usefull if we only // want to compute the viscosity solution of the considered equation by // the Fast Marching Method. Nevertheless optimal Dymamics is usefull in // many application as for example fibers tracking in DTI (see papers // associated to this code). setOptDynamics(x,y,z, optDymamics1[gs1][gs2][gs3],optDymamics2[gs1][gs2][gs3],optDymamics3[gs1][gs2][gs3]); // we can then return the value for update... return minRoot; }// End of the methode "_UpdateValue". }; // End of the definition of the class.} // End of namespace FastLevelSet.#endif?? 快捷鍵說明
復(fù)制代碼
Ctrl + C
搜索代碼
Ctrl + F
全屏模式
F11
切換主題
Ctrl + Shift + D
顯示快捷鍵
?
增大字號
Ctrl + =
減小字號
Ctrl + -