diff --git a/src/include/tra_subs.h b/src/include/tra_subs.h new file mode 100644 index 0000000000000000000000000000000000000000..f4f50ecd690474032a7e2e5795a9ef5349f9cfb8 --- /dev/null +++ b/src/include/tra_subs.h @@ -0,0 +1,65 @@ +#include <fstream> +#include <complex> + +#ifndef INCLUDE_TRA_SUBS_H_ +#define INCLUDE_TRA_SUBS_H_ +#endif + +/*! C++ porting of FRFMER + * + * \param nkv: `int` QUESTION: definition? + * \param vkm: `double` QUESTION: definition? + * \param vkv: `double *` QUESTION: definition? + * \param vknmx: `double` QUESTION: definition? + * \param apfafa: `double` QUESTION: definition? + * \param tra: `double` QUESTION: definition? + * \param spd: `double` QUESTION: definition? + * \param rir: `double` QUESTION: definition? + * \param fctn: `double` QUESTION: definition? + * \param le: `int` QUESTION: definition? + * \param lmode: `int` QUESTION: definition? + * \param pmf: `double` QUESTION: definition? + * \param tt1: `fstream &` Handle to first temporary binary file. + * \param tt2: `fstream &` Handle to second temporary binary file. + */ +void frfmer( + int nkv, double vkm, double *vkv, double vknmx, double apfafa, double tra, + double spd, double rir, double ftcn, int le, int lmode, double pmf, + std::fstream &tt1, std::fstream &tt2 +) { + const int nlemt = le * (le + 2) * 2; + const std::complex<double> cc0(0.0, 0.0); + std::complex<double> *wk = new std::complex<double>[nlemt](); + double sq = vkm * vkm; + for (int jy90 = 0; jy90 < nkv; jy90++) { + double vky = vkv[jy90]; + double sqy = vky * vky; + for (int jx80 = 0; jx80 < nkv; jx80++) { + double vkx = vkv[jx80]; + double sqx = vkx * vkx; + double sqn = sqx + sqy; + double vkn = sqrt(sqn); + if (vkn <= vknmx) { + double vkz = sqrt(sq - sqn); + // Would call WAMFF(WK,VKX,VKY,VKZ,LE,APFAFA,TRA,SPD,RIR,FTCN,LMODE,PMF) + for (int j = 0; j < nlemt; j++) { + double vreal = wk[j].real(); + double vimag = wk[j].imag(); + tt1.write(reinterpret_cast<char *>(&vreal), sizeof(double)); + tt1.write(reinterpret_cast<char *>(&vimag), sizeof(double)); + } + tt2.write(reinterpret_cast<char *>(&vkz), sizeof(double)); + } else { // label 50 + for (int j = 0; j < nlemt; j++) { + double vreal = 0.0; + double vimag = 0.0; + tt1.write(reinterpret_cast<char *>(&vreal), sizeof(double)); + tt1.write(reinterpret_cast<char *>(&vimag), sizeof(double)); + } + double vkz = 0.0; + tt2.write(reinterpret_cast<char *>(&vkz), sizeof(double)); + } + } // jx80 loop + } // jy90 loop + delete[] wk; +} diff --git a/src/trapping/frfme.cpp b/src/trapping/frfme.cpp index b34cb37ba6147deb187b2035fc155ba01cd62a2a..5eaf03eeb109f09066773e568b400efa22c1f604 100644 --- a/src/trapping/frfme.cpp +++ b/src/trapping/frfme.cpp @@ -9,6 +9,9 @@ //#ifndef INCLUDE_SPH_SUBS_H_ //#include "../include/sph_subs.h" //#endif +#ifndef INCLUDE_TRA_SUBS_H_ +#include "../include/tra_subs.h" +#endif using namespace std; @@ -293,7 +296,7 @@ int main() { temptape2.open(temp_name2.c_str(), ios::out | ios::binary); for (int jx = 0; jx < nkv; jx++) temptape2.write(reinterpret_cast<char *>(&(vkv[jx])), sizeof(double)); - // Would call FRFMER(NKV,VKM,VKV,VKNMX,APFAFA,TRA,SPD,RIR,FTCN,LM,LMODE,PMF,ITT1,ITT2) + frfmer(nkv, vkm, vkv, vknmx, apfafa, tra, spd, rir, ftcn, lm, lmode, pmf, temptape1, temptape2); temptape1.write(reinterpret_cast<char *>(&apfafa), sizeof(double)); // Place-holder temptape1.close(); temptape2.write(reinterpret_cast<char *>(&apfafa), sizeof(double)); @@ -339,5 +342,6 @@ int main() { for (int wsi = wsum_size - 1; wsi > -1; wsi--) delete[] wsum[wsi]; delete[] wsum; } + printf("Done.\n"); return 0; }