Loading src/trapping/cfrfme.cpp +83 −129 Original line number Original line Diff line number Diff line Loading @@ -4,6 +4,7 @@ */ */ #include <complex> #include <complex> #include <cstdio> #include <cstdio> #include <exception> #include <fstream> #include <fstream> #include <regex> #include <regex> #include <string> #include <string> Loading @@ -16,10 +17,18 @@ #include "../include/Commons.h" #include "../include/Commons.h" #endif #endif #ifndef INCLUDE_CONFIGURATION_H_ #include "../include/Configuration.h" #endif #ifndef INCLUDE_SPH_SUBS_H_ #ifndef INCLUDE_SPH_SUBS_H_ #include "../include/sph_subs.h" #include "../include/sph_subs.h" #endif #endif #ifndef INCLUDE_TFRFME_H_ #include "../include/tfrfme.h" #endif #ifndef INCLUDE_TRA_SUBS_H_ #ifndef INCLUDE_TRA_SUBS_H_ #include "../include/tra_subs.h" #include "../include/tra_subs.h" #endif #endif Loading @@ -32,13 +41,13 @@ using namespace std; * \param output_path: `string` Directory to write the output files in. * \param output_path: `string` Directory to write the output files in. */ */ void frfme(string data_file, string output_path) { void frfme(string data_file, string output_path) { string tfrfme_name = output_path + "/c_TFRFME"; string tfrfme_name = output_path + "/c_TFRFME.hd5"; fstream tfrfme; //fstream tfrfme; TFRFME *tfrfme = NULL; char namef[7]; char namef[7]; char more; char more; double *xv = NULL, *yv = NULL, *zv = NULL; double *vkv = NULL, **vkzm = NULL; double *vkv = NULL, **vkzm = NULL; complex<double> *wk = NULL, **w = NULL, **wsum = NULL; complex<double> *wk = NULL, **w = NULL; const complex<double> cc0(0.0, 0.0); const complex<double> cc0(0.0, 0.0); const complex<double> uim(0.0, 1.0); const complex<double> uim(0.0, 1.0); int line_count = 0, last_read_line = 0; int line_count = 0, last_read_line = 0; Loading @@ -64,29 +73,23 @@ void frfme(string data_file, string output_path) { // End of vector size variables // End of vector size variables if (jlmf != 1) { if (jlmf != 1) { int nxv, nyv, nzv; int nxv, nyv, nzv; tfrfme.open(tfrfme_name, ios::in | ios::binary); if (tfrfme == NULL) tfrfme = TFRFME::from_binary(tfrfme_name, "HDF5"); if (tfrfme.is_open()) { if (tfrfme != NULL) { tfrfme.read(reinterpret_cast<char *>(&lmode), sizeof(int)); lmode = (int)tfrfme->get_param("lmode"); tfrfme.read(reinterpret_cast<char *>(&lm), sizeof(int)); lm = (int)tfrfme->get_param("lm"); tfrfme.read(reinterpret_cast<char *>(&nkv), sizeof(int)); nkv = (int)tfrfme->get_param("nkv"); tfrfme.read(reinterpret_cast<char *>(&nxv), sizeof(int)); nxv = (int)tfrfme->get_param("nkv"); tfrfme.read(reinterpret_cast<char *>(&nyv), sizeof(int)); nyv = (int)tfrfme->get_param("nkv"); tfrfme.read(reinterpret_cast<char *>(&nzv), sizeof(int)); nzv = (int)tfrfme->get_param("nkv"); tfrfme.read(reinterpret_cast<char *>(&vk), sizeof(double)); vk = tfrfme->get_param("vk"); tfrfme.read(reinterpret_cast<char *>(&exri), sizeof(double)); exri = tfrfme->get_param("exri"); tfrfme.read(reinterpret_cast<char *>(&an), sizeof(double)); an = tfrfme->get_param("an"); tfrfme.read(reinterpret_cast<char *>(&ff), sizeof(double)); ff = tfrfme->get_param("ff"); tfrfme.read(reinterpret_cast<char *>(&tra), sizeof(double)); tra = tfrfme->get_param("tra"); tfrfme.read(reinterpret_cast<char *>(&spd), sizeof(double)); spd = tfrfme->get_param("spd"); tfrfme.read(reinterpret_cast<char *>(&frsh), sizeof(double)); frsh = tfrfme->get_param("frsh"); tfrfme.read(reinterpret_cast<char *>(&exril), sizeof(double)); exril = tfrfme->get_param("exril"); vkv = new double[nkv](); vkv = new double[nkv](); xv = new double[nxv](); yv = new double[nyv](); zv = new double[nzv](); for (int xi = 0; xi < nxv; xi++) tfrfme.read(reinterpret_cast<char *>(&(xv[xi])), sizeof(double)); for (int yi = 0; yi < nxv; yi++) tfrfme.read(reinterpret_cast<char *>(&(yv[yi])), sizeof(double)); for (int zi = 0; zi < nxv; zi++) tfrfme.read(reinterpret_cast<char *>(&(zv[zi])), sizeof(double)); fstream temptape2; fstream temptape2; string tempname2 = output_path + "/c_TEMPTAPE2"; string tempname2 = output_path + "/c_TEMPTAPE2"; temptape2.open(tempname2.c_str(), ios::in | ios::binary); temptape2.open(tempname2.c_str(), ios::in | ios::binary); Loading Loading @@ -116,15 +119,6 @@ void frfme(string data_file, string output_path) { } else { } else { printf("ERROR: could not open TEMPTAPE2 file.\n"); printf("ERROR: could not open TEMPTAPE2 file.\n"); } } for (int ixyz12 = 0; ixyz12 < nrvc; ixyz12++) { for (int j12 = 0; j12 < jlmf - 1; j12++) { double vreal, vimag; tfrfme.read(reinterpret_cast<char *>(&vreal), sizeof(double)); tfrfme.read(reinterpret_cast<char *>(&vimag), sizeof(double)); wsum[j12][ixyz12] = complex<double>(vreal, vimag); } // j12 loop } // ixyz12 loop tfrfme.close(); } else { } else { printf("ERROR: could not open TFRFME file.\n"); printf("ERROR: could not open TFRFME file.\n"); } } Loading Loading @@ -179,23 +173,22 @@ void frfme(string data_file, string output_path) { re = regex("[0-9]+"); re = regex("[0-9]+"); regex_search(str_target, m, re); regex_search(str_target, m, re); int ixi = stoi(m.str()); int ixi = stoi(m.str()); fstream tedf; string tedf_name = output_path + "/" + namef + ".hd5"; string tedf_name = output_path + "/" + namef; ScattererConfiguration *tedf = ScattererConfiguration::from_binary(tedf_name, "HDF5"); tedf.open(tedf_name.c_str(), ios::in | ios::binary); if (tedf != NULL) { if (tedf.is_open()) { int iduml, idum; int iduml, idum; tedf.read(reinterpret_cast<char *>(&iduml), sizeof(int)); iduml = (int)tedf->get_param("nsph"); for (int i = 0; i < iduml; i++) tedf.read(reinterpret_cast<char *>(&idum), sizeof(int)); idum = tedf->get_iog(iduml - 1); tedf.read(reinterpret_cast<char *>(&exdc), sizeof(double)); exdc = tedf->get_param("exdc"); tedf.read(reinterpret_cast<char *>(&wp), sizeof(double)); wp = tedf->get_param("wp"); tedf.read(reinterpret_cast<char *>(&xip), sizeof(double)); xip = tedf->get_param("xip"); tedf.read(reinterpret_cast<char *>(&idfc), sizeof(int)); idfc = (int)tedf->get_param("idfc"); tedf.read(reinterpret_cast<char *>(&nxi), sizeof(int)); nxi = (int)tedf->get_param("nxi"); if (idfc >= 0) { if (idfc >= 0) { if (ixi <= nxi) { if (ixi <= nxi) { for (int i = 0; i < ixi; i++) tedf.read(reinterpret_cast<char *>(&xi), sizeof(double)); xi = tedf->get_scale(ixi - 1); } else { // label 96 } else { // label 96 tedf.close(); delete tedf; // label 98 // label 98 string output_name = output_path + "/c_OFRFME"; string output_name = output_path + "/c_OFRFME"; FILE *output = fopen(output_name.c_str(), "w"); FILE *output = fopen(output_name.c_str(), "w"); Loading @@ -206,7 +199,7 @@ void frfme(string data_file, string output_path) { xi = xip; xi = xip; } } // label 20 // label 20 tedf.close(); delete tedf; double wn = wp / 3.0e8; double wn = wp / 3.0e8; vk = xi * wn; vk = xi * wn; exri = sqrt(exdc); exri = sqrt(exdc); Loading Loading @@ -244,26 +237,24 @@ void frfme(string data_file, string output_path) { int nxs = nxsh * 2; int nxs = nxsh * 2; int nxv = nxs + 1; int nxv = nxs + 1; int nxshpo = nxsh + 1; int nxshpo = nxsh + 1; xv = new double[nxv](); for (int i24 = nxshpo; i24 <= nxs; i24++) { xv[i24] = xv[i24 - 1] + delxyz; xv[nxv - i24 - 1] = -xv[i24]; } // i24 loop int nys = nysh * 2; int nys = nysh * 2; int nyv = nys + 1; int nyv = nys + 1; int nyshpo = nysh + 1; int nyshpo = nysh + 1; yv = new double[nyv](); for (int i25 = nyshpo; i25 <= nys; i25++) { yv[i25] = yv[i25 - 1] + delxyz; yv[nyv - i25 - 1] = -yv[i25]; } // i25 loop int nzs = nzsh * 2; int nzs = nzsh * 2; int nzv = nzs + 1; int nzv = nzs + 1; int nzshpo = nzsh + 1; int nzshpo = nzsh + 1; zv = new double[nzv](); tfrfme = new TFRFME(lmode, lm, nkv, nxv, nyv, nzv); for (int i24 = nxshpo; i24 <= nxs; i24++) { tfrfme->set_x(i24, tfrfme->get_x(i24 - 1) + delxyz); tfrfme->set_x(nxv - i24 - 1, -tfrfme->get_x(i24)); } // i24 loop for (int i25 = nyshpo; i25 <= nys; i25++) { tfrfme->set_y(i25, tfrfme->get_y(i25 - 1) + delxyz); tfrfme->set_y(nyv - i25 - 1, -tfrfme->get_y(i25)); } // i25 loop for (int i27 = nzshpo; i27 <= nzs; i27++) { for (int i27 = nzshpo; i27 <= nzs; i27++) { zv[i27] = zv[i27 - 1] + delxyz; tfrfme->set_z(i27, tfrfme->get_z(i27 - 1) + delxyz); zv[nzv - i27 - 1] = -zv[i27]; tfrfme->set_z(nzv - i27 - 1, -tfrfme->get_z(i27)); } // i27 loop } // i27 loop int nrvc = nxv * nyv * nzv; int nrvc = nxv * nyv * nzv; int nkshpo = nksh + 1; int nkshpo = nksh + 1; Loading @@ -271,28 +262,15 @@ void frfme(string data_file, string output_path) { vkv[i28] = vkv[i28 - 1] + delk; vkv[i28] = vkv[i28 - 1] + delk; vkv[nkv - i28 - 1] = -vkv[i28]; vkv[nkv - i28 - 1] = -vkv[i28]; } // i28 loop } // i28 loop tfrfme.open(tfrfme_name.c_str(), ios::out | ios::binary); if (tfrfme != NULL) { if (tfrfme.is_open()) { tfrfme->set_param("vk", vk); tfrfme.write(reinterpret_cast<char *>(&lmode), sizeof(int)); tfrfme->set_param("exri", exri); tfrfme.write(reinterpret_cast<char *>(&lm), sizeof(int)); tfrfme->set_param("an", an); tfrfme.write(reinterpret_cast<char *>(&nkv), sizeof(int)); tfrfme->set_param("ff", ff); tfrfme.write(reinterpret_cast<char *>(&nxv), sizeof(int)); tfrfme->set_param("tra", tra); tfrfme.write(reinterpret_cast<char *>(&nyv), sizeof(int)); tfrfme->set_param("spd", spd); tfrfme.write(reinterpret_cast<char *>(&nzv), sizeof(int)); tfrfme->set_param("frsh", frsh); tfrfme.write(reinterpret_cast<char *>(&vk), sizeof(double)); tfrfme->set_param("exril", exril); tfrfme.write(reinterpret_cast<char *>(&exri), sizeof(double)); tfrfme.write(reinterpret_cast<char *>(&an), sizeof(double)); tfrfme.write(reinterpret_cast<char *>(&ff), sizeof(double)); tfrfme.write(reinterpret_cast<char *>(&tra), sizeof(double)); tfrfme.write(reinterpret_cast<char *>(&spd), sizeof(double)); tfrfme.write(reinterpret_cast<char *>(&frsh), sizeof(double)); tfrfme.write(reinterpret_cast<char *>(&exril), sizeof(double)); for (int xi = 0; xi < nxv; xi++) tfrfme.write(reinterpret_cast<char *>(&(xv[xi])), sizeof(double)); for (int yi = 0; yi < nyv; yi++) tfrfme.write(reinterpret_cast<char *>(&(yv[yi])), sizeof(double)); for (int zi = 0; zi < nzv; zi++) tfrfme.write(reinterpret_cast<char *>(&(zv[zi])), sizeof(double)); fstream temptape1, temptape2; fstream temptape1, temptape2; string temp_name1 = output_path + "/c_TEMPTAPE1"; string temp_name1 = output_path + "/c_TEMPTAPE1"; string temp_name2 = output_path + "/c_TEMPTAPE2"; string temp_name2 = output_path + "/c_TEMPTAPE2"; Loading Loading @@ -330,11 +308,6 @@ void frfme(string data_file, string output_path) { } } } // jy40 loop } // jy40 loop temptape2.close(); temptape2.close(); if (wsum != NULL) { for (int wsi = wsum_size - 1; wsi > -1; wsi--) delete[] wsum[wsi]; delete[] wsum; } wsum = new complex<double>*[nlmmt]; for (int j80 = jlmf - 1; j80 < jlml; j80++) { for (int j80 = jlmf - 1; j80 < jlml; j80++) { // w matrix // w matrix if (w != NULL) { if (w != NULL) { Loading @@ -359,13 +332,13 @@ void frfme(string data_file, string output_path) { } // jy50 loop } // jy50 loop temptape1.close(); temptape1.close(); int ixyz = 0; int ixyz = 0; wsum[j80] = new complex<double>[nrvc](); for (int wj = 0; wj < nrvc; wj++) tfrfme->set_matrix_element(j80, wj, cc0); for (int iz75 = 0; iz75 < nzv; iz75++) { for (int iz75 = 0; iz75 < nzv; iz75++) { double z = zv[iz75] + frsh; double z = tfrfme->get_z(iz75) + frsh; for (int iy70 = 0; iy70 < nyv; iy70++) { for (int iy70 = 0; iy70 < nyv; iy70++) { double y = yv[iy70]; double y = tfrfme->get_y(iy70); for (int ix65 = 0; ix65 < nxv; ix65++) { for (int ix65 = 0; ix65 < nxv; ix65++) { double x = xv[ix65]; double x = tfrfme->get_x(ix65); ixyz++; ixyz++; complex<double> sumy = cc0; complex<double> sumy = cc0; for (int jy60 = 0; jy60 < nkv; jy60++) { for (int jy60 = 0; jy60 < nkv; jy60++) { Loading @@ -385,41 +358,29 @@ void frfme(string data_file, string output_path) { if (jy60 == 0 || jy60 == nkv - 1) sumx *= 0.5; if (jy60 == 0 || jy60 == nkv - 1) sumx *= 0.5; sumy += sumx; sumy += sumx; } // jy60 loop } // jy60 loop wsum[j80][ixyz - 1] = sumy * delks; tfrfme->set_matrix_element(j80, ixyz - 1, sumy * delks); } // ix65 loop } // ix65 loop } // iy70 loop } // iy70 loop } // iz75 loop } // iz75 loop } // j80 loop } // j80 loop if (jlmf != 1) { if (jlmf != 1) { tfrfme.open(tfrfme_name, ios::in | ios::out | ios::binary); lmode = (int)tfrfme->get_param("lmode"); tfrfme.read(reinterpret_cast<char *>(&lmode), sizeof(int)); lm = (int)tfrfme->get_param("lm"); tfrfme.read(reinterpret_cast<char *>(&lm), sizeof(int)); nkv = (int)tfrfme->get_param("nkv"); tfrfme.read(reinterpret_cast<char *>(&nkv), sizeof(int)); nxv = (int)tfrfme->get_param("nxv"); tfrfme.read(reinterpret_cast<char *>(&nxv), sizeof(int)); nyv = (int)tfrfme->get_param("nyv"); tfrfme.read(reinterpret_cast<char *>(&nyv), sizeof(int)); nzv = (int)tfrfme->get_param("nzv"); tfrfme.read(reinterpret_cast<char *>(&nzv), sizeof(int)); vk = tfrfme->get_param("vk"); tfrfme.read(reinterpret_cast<char *>(&vk), sizeof(double)); exri = tfrfme->get_param("exri"); tfrfme.read(reinterpret_cast<char *>(&exri), sizeof(double)); an = tfrfme->get_param("an"); tfrfme.read(reinterpret_cast<char *>(&an), sizeof(double)); ff = tfrfme->get_param("ff"); tfrfme.read(reinterpret_cast<char *>(&ff), sizeof(double)); tra = tfrfme->get_param("tra"); tfrfme.read(reinterpret_cast<char *>(&tra), sizeof(double)); spd = tfrfme->get_param("spd"); tfrfme.read(reinterpret_cast<char *>(&spd), sizeof(double)); frsh = tfrfme->get_param("frsh"); tfrfme.read(reinterpret_cast<char *>(&frsh), sizeof(double)); exril = tfrfme->get_param("exril"); tfrfme.read(reinterpret_cast<char *>(&exril), sizeof(double)); for (int i = 0; i < nxv; i++) tfrfme.read(reinterpret_cast<char *>(&(xv[i])), sizeof(double)); for (int i = 0; i < nyv; i++) tfrfme.read(reinterpret_cast<char *>(&(yv[i])), sizeof(double)); for (int i = 0; i < nzv; i++) tfrfme.read(reinterpret_cast<char *>(&(zv[i])), sizeof(double)); } } // label 88 // label 88 for (int ixyz = 0; ixyz < nrvc; ixyz++) { tfrfme->write_binary(tfrfme_name, "HDF5"); for (int j = 0; j < jlml; j++) { double vreal = wsum[j][ixyz].real(); double vimag = wsum[j][ixyz].imag(); tfrfme.write(reinterpret_cast<char *>(&vreal), sizeof(double)); tfrfme.write(reinterpret_cast<char *>(&vimag), sizeof(double)); } // j loop } // ixyz loop tfrfme.close(); string output_name = output_path + "/c_OFRFME"; string output_name = output_path + "/c_OFRFME"; FILE *output = fopen(output_name.c_str(), "w"); FILE *output = fopen(output_name.c_str(), "w"); fprintf(output, " IF JLML < NLMMT, PRESERVE TEMPTAPE1, TEMPTAPE2, AND TFRFRME,\n"); fprintf(output, " IF JLML < NLMMT, PRESERVE TEMPTAPE1, TEMPTAPE2, AND TFRFRME,\n"); Loading @@ -441,11 +402,8 @@ void frfme(string data_file, string output_path) { } } } } // label 45 // label 45 if (tfrfme.is_open()) tfrfme.close(); if (tfrfme != NULL) delete tfrfme; delete[] file_lines; delete[] file_lines; if (xv != NULL) delete[] xv; if (yv != NULL) delete[] yv; if (zv != NULL) delete[] zv; if (vkv != NULL) delete[] vkv; if (vkv != NULL) delete[] vkv; if (vkzm != NULL) { if (vkzm != NULL) { for (int vki = nkv - 1; vki > -1; vki--) delete[] vkzm[vki]; for (int vki = nkv - 1; vki > -1; vki--) delete[] vkzm[vki]; Loading @@ -455,10 +413,6 @@ void frfme(string data_file, string output_path) { for (int wi = nkv - 1; wi > -1; wi--) delete[] w[wi]; for (int wi = nkv - 1; wi > -1; wi--) delete[] w[wi]; delete[] w; delete[] w; } } if (wsum != NULL) { for (int wsi = wsum_size - 1; wsi > -1; wsi--) delete[] wsum[wsi]; delete[] wsum; } if (wk != NULL) delete[] wk; if (wk != NULL) delete[] wk; printf("Done.\n"); printf("Done.\n"); } } src/trapping/clffft.cpp +58 −46 Original line number Original line Diff line number Diff line Loading @@ -4,6 +4,7 @@ */ */ #include <complex> #include <complex> #include <cstdio> #include <cstdio> #include <exception> #include <fstream> #include <fstream> #include <regex> #include <regex> #include <string> #include <string> Loading @@ -20,6 +21,10 @@ #include "../include/sph_subs.h" #include "../include/sph_subs.h" #endif #endif #ifndef INCLUDE_TFRFME_H_ #include "../include/tfrfme.h" #endif #ifndef INCLUDE_TRA_SUBS_H_ #ifndef INCLUDE_TRA_SUBS_H_ #include "../include/tra_subs.h" #include "../include/tra_subs.h" #endif #endif Loading @@ -38,7 +43,6 @@ void lffft(string data_file, string output_path) { fstream tlfff, tlfft; fstream tlfff, tlfft; double ****zpv = NULL; double ****zpv = NULL; double *xv = NULL, *yv = NULL, *zv = NULL; complex<double> *ac = NULL, *ws = NULL, *wsl = NULL; complex<double> *ac = NULL, *ws = NULL, *wsl = NULL; complex<double> **am0m = NULL; complex<double> **am0m = NULL; complex<double> **amd = NULL; complex<double> **amd = NULL; Loading Loading @@ -140,37 +144,31 @@ void lffft(string data_file, string output_path) { } } // label 150 // label 150 ttms.close(); ttms.close(); fstream binary_input; TFRFME *tfrfme = NULL; string binary_name; string binary_name; if (jss != 1) binary_name = output_path + "/c_TFRFME"; if (jss != 1) binary_name = output_path + "/c_TFRFME.hd5"; else binary_name = output_path + "/c_TWS"; else binary_name = output_path + "/c_TWS"; binary_input.open(binary_name, ios::in | ios::binary); tfrfme = TFRFME::from_binary(binary_name, "HDF5"); if (binary_input.is_open()) { if (tfrfme != NULL) { int lmode, lm, nkv, nxv, nyv, nzv; int lmode, lm, nkv, nxv, nyv, nzv; double vk, exri, an, ff, tra; double vk, exri, an, ff, tra; double spd, frsh, exril; double spd, frsh, exril; binary_input.read(reinterpret_cast<char *>(&lmode), sizeof(int)); lmode = (int)tfrfme->get_param("lmode"); binary_input.read(reinterpret_cast<char *>(&lm), sizeof(int)); lm = (int)tfrfme->get_param("lm"); binary_input.read(reinterpret_cast<char *>(&nkv), sizeof(int)); nkv = (int)tfrfme->get_param("nkv"); binary_input.read(reinterpret_cast<char *>(&nxv), sizeof(int)); nxv = (int)tfrfme->get_param("nxv"); binary_input.read(reinterpret_cast<char *>(&nyv), sizeof(int)); nyv = (int)tfrfme->get_param("nyv"); binary_input.read(reinterpret_cast<char *>(&nzv), sizeof(int)); nzv = (int)tfrfme->get_param("nzv"); if (lm >= le) { if (lm >= le) { binary_input.read(reinterpret_cast<char *>(&vk), sizeof(double)); vk = tfrfme->get_param("vk"); binary_input.read(reinterpret_cast<char *>(&exri), sizeof(double)); exri = tfrfme->get_param("exri"); binary_input.read(reinterpret_cast<char *>(&an), sizeof(double)); an = tfrfme->get_param("an"); binary_input.read(reinterpret_cast<char *>(&ff), sizeof(double)); ff = tfrfme->get_param("ff"); binary_input.read(reinterpret_cast<char *>(&tra), sizeof(double)); tra = tfrfme->get_param("tra"); if (vk == vks && exri == exris) { if (vk == vks && exri == exris) { binary_input.read(reinterpret_cast<char *>(&spd), sizeof(double)); spd = tfrfme->get_param("spd"); binary_input.read(reinterpret_cast<char *>(&frsh), sizeof(double)); frsh = tfrfme->get_param("frsh"); binary_input.read(reinterpret_cast<char *>(&exril), sizeof(double)); exril = tfrfme->get_param("exril"); xv = new double[nxv]; for (int i = 0; i < nxv; i++) binary_input.read(reinterpret_cast<char *>(&(xv[i])), sizeof(double)); yv = new double[nyv]; for (int i = 0; i < nyv; i++) binary_input.read(reinterpret_cast<char *>(&(yv[i])), sizeof(double)); zv = new double[nzv]; for (int i = 0; i < nzv; i++) binary_input.read(reinterpret_cast<char *>(&(zv[i])), sizeof(double)); bool goto160 = false; bool goto160 = false; if (jft <= 0) { if (jft <= 0) { zpv = new double***[le]; zpv = new double***[le]; Loading Loading @@ -204,12 +202,18 @@ void lffft(string data_file, string output_path) { tlfff.write(reinterpret_cast<char *>(&spd), sizeof(double)); tlfff.write(reinterpret_cast<char *>(&spd), sizeof(double)); tlfff.write(reinterpret_cast<char *>(&frsh), sizeof(double)); tlfff.write(reinterpret_cast<char *>(&frsh), sizeof(double)); tlfff.write(reinterpret_cast<char *>(&exril), sizeof(double)); tlfff.write(reinterpret_cast<char *>(&exril), sizeof(double)); for (int i = 0; i < nxv; i++) for (int i = 0; i < nxv; i++) { tlfff.write(reinterpret_cast<char *>(&(xv[i])), sizeof(double)); double x = tfrfme->get_x(i); for (int i = 0; i < nyv; i++) tlfff.write(reinterpret_cast<char *>(&x), sizeof(double)); tlfff.write(reinterpret_cast<char *>(&(yv[i])), sizeof(double)); } for (int i = 0; i < nzv; i++) for (int i = 0; i < nyv; i++) { tlfff.write(reinterpret_cast<char *>(&(zv[i])), sizeof(double)); double y = tfrfme->get_y(i); tlfff.write(reinterpret_cast<char *>(&y), sizeof(double)); } for (int i = 0; i < nzv; i++) { double z = tfrfme->get_z(i); tlfff.write(reinterpret_cast<char *>(&z), sizeof(double)); } if (jft < 0) goto160 = true; if (jft < 0) goto160 = true; } else { // Should never happen. } else { // Should never happen. printf("ERROR: could not open TLFFF file.\n"); printf("ERROR: could not open TLFFF file.\n"); Loading @@ -236,12 +240,18 @@ void lffft(string data_file, string output_path) { tlfft.write(reinterpret_cast<char *>(&spd), sizeof(double)); tlfft.write(reinterpret_cast<char *>(&spd), sizeof(double)); tlfft.write(reinterpret_cast<char *>(&frsh), sizeof(double)); tlfft.write(reinterpret_cast<char *>(&frsh), sizeof(double)); tlfft.write(reinterpret_cast<char *>(&exril), sizeof(double)); tlfft.write(reinterpret_cast<char *>(&exril), sizeof(double)); for (int i = 0; i < nxv; i++) for (int i = 0; i < nxv; i++) { tlfft.write(reinterpret_cast<char *>(&(xv[i])), sizeof(double)); double x = tfrfme->get_x(i); for (int i = 0; i < nyv; i++) tlfft.write(reinterpret_cast<char *>(&x), sizeof(double)); tlfft.write(reinterpret_cast<char *>(&(yv[i])), sizeof(double)); } for (int i = 0; i < nzv; i++) for (int i = 0; i < nyv; i++) { tlfft.write(reinterpret_cast<char *>(&(zv[i])), sizeof(double)); double y = tfrfme->get_y(i); tlfft.write(reinterpret_cast<char *>(&y), sizeof(double)); } for (int i = 0; i < nzv; i++) { double z = tfrfme->get_z(i); tlfft.write(reinterpret_cast<char *>(&z), sizeof(double)); } } else { // Should never happen. } else { // Should never happen. printf("ERROR: could not open TLFFT file.\n"); printf("ERROR: could not open TLFFT file.\n"); } } Loading @@ -262,13 +272,18 @@ void lffft(string data_file, string output_path) { for (int iy475 = 0; iy475 < nyv; iy475++) { for (int iy475 = 0; iy475 < nyv; iy475++) { for (int ix475 = 0; ix475 < nxv; ix475++) { for (int ix475 = 0; ix475 < nxv; ix475++) { for (int i = 0; i < nlmmt; i++) { for (int i = 0; i < nlmmt; i++) { double vreal, vimag; //double vreal, vimag; binary_input.read(reinterpret_cast<char *>(&vreal), sizeof(double)); //binary_input.read(reinterpret_cast<char *>(&vreal), sizeof(double)); binary_input.read(reinterpret_cast<char *>(&vimag), sizeof(double)); //binary_input.read(reinterpret_cast<char *>(&vimag), sizeof(double)); int row = i; int col = (nzv * iz475) + (nyv * iy475) + ix475; complex<double> value = tfrfme->get_matrix_element(row, col); if (lm <= le) { if (lm <= le) { ws[i] = complex<double>(vreal, vimag); //ws[i] = complex<double>(vreal, vimag); ws[i] = value; } else { // label 170 } else { // label 170 wsl[i] = complex<double>(vreal, vimag); //wsl[i] = complex<double>(vreal, vimag); wsl[i] = value; for (int i175 = 0; i175 < cil->nlem; i175++) { for (int i175 = 0; i175 < cil->nlem; i175++) { int ie = i175 + cil->nlem; int ie = i175 + cil->nlem; int iel = i175 + nlmm; int iel = i175 + nlmm; Loading Loading @@ -384,7 +399,7 @@ void lffft(string data_file, string output_path) { fclose(output67); fclose(output67); } } } } binary_input.close(); delete tfrfme; } else { } else { printf("ERROR: could not open binary input file %s.\n", binary_name.c_str()); printf("ERROR: could not open binary input file %s.\n", binary_name.c_str()); } } Loading @@ -395,9 +410,6 @@ void lffft(string data_file, string output_path) { // Clean up memory // Clean up memory if (ac != NULL) delete[] ac; if (ac != NULL) delete[] ac; if (ws != NULL) delete[] ws; if (ws != NULL) delete[] ws; if (xv != NULL) delete[] xv; if (yv != NULL) delete[] yv; if (zv != NULL) delete[] zv; if (wsl != NULL) delete[] wsl; if (wsl != NULL) delete[] wsl; if (tmsm != NULL) delete[] tmsm; if (tmsm != NULL) delete[] tmsm; if (tmse != NULL) delete[] tmse; if (tmse != NULL) delete[] tmse; Loading Loading
src/trapping/cfrfme.cpp +83 −129 Original line number Original line Diff line number Diff line Loading @@ -4,6 +4,7 @@ */ */ #include <complex> #include <complex> #include <cstdio> #include <cstdio> #include <exception> #include <fstream> #include <fstream> #include <regex> #include <regex> #include <string> #include <string> Loading @@ -16,10 +17,18 @@ #include "../include/Commons.h" #include "../include/Commons.h" #endif #endif #ifndef INCLUDE_CONFIGURATION_H_ #include "../include/Configuration.h" #endif #ifndef INCLUDE_SPH_SUBS_H_ #ifndef INCLUDE_SPH_SUBS_H_ #include "../include/sph_subs.h" #include "../include/sph_subs.h" #endif #endif #ifndef INCLUDE_TFRFME_H_ #include "../include/tfrfme.h" #endif #ifndef INCLUDE_TRA_SUBS_H_ #ifndef INCLUDE_TRA_SUBS_H_ #include "../include/tra_subs.h" #include "../include/tra_subs.h" #endif #endif Loading @@ -32,13 +41,13 @@ using namespace std; * \param output_path: `string` Directory to write the output files in. * \param output_path: `string` Directory to write the output files in. */ */ void frfme(string data_file, string output_path) { void frfme(string data_file, string output_path) { string tfrfme_name = output_path + "/c_TFRFME"; string tfrfme_name = output_path + "/c_TFRFME.hd5"; fstream tfrfme; //fstream tfrfme; TFRFME *tfrfme = NULL; char namef[7]; char namef[7]; char more; char more; double *xv = NULL, *yv = NULL, *zv = NULL; double *vkv = NULL, **vkzm = NULL; double *vkv = NULL, **vkzm = NULL; complex<double> *wk = NULL, **w = NULL, **wsum = NULL; complex<double> *wk = NULL, **w = NULL; const complex<double> cc0(0.0, 0.0); const complex<double> cc0(0.0, 0.0); const complex<double> uim(0.0, 1.0); const complex<double> uim(0.0, 1.0); int line_count = 0, last_read_line = 0; int line_count = 0, last_read_line = 0; Loading @@ -64,29 +73,23 @@ void frfme(string data_file, string output_path) { // End of vector size variables // End of vector size variables if (jlmf != 1) { if (jlmf != 1) { int nxv, nyv, nzv; int nxv, nyv, nzv; tfrfme.open(tfrfme_name, ios::in | ios::binary); if (tfrfme == NULL) tfrfme = TFRFME::from_binary(tfrfme_name, "HDF5"); if (tfrfme.is_open()) { if (tfrfme != NULL) { tfrfme.read(reinterpret_cast<char *>(&lmode), sizeof(int)); lmode = (int)tfrfme->get_param("lmode"); tfrfme.read(reinterpret_cast<char *>(&lm), sizeof(int)); lm = (int)tfrfme->get_param("lm"); tfrfme.read(reinterpret_cast<char *>(&nkv), sizeof(int)); nkv = (int)tfrfme->get_param("nkv"); tfrfme.read(reinterpret_cast<char *>(&nxv), sizeof(int)); nxv = (int)tfrfme->get_param("nkv"); tfrfme.read(reinterpret_cast<char *>(&nyv), sizeof(int)); nyv = (int)tfrfme->get_param("nkv"); tfrfme.read(reinterpret_cast<char *>(&nzv), sizeof(int)); nzv = (int)tfrfme->get_param("nkv"); tfrfme.read(reinterpret_cast<char *>(&vk), sizeof(double)); vk = tfrfme->get_param("vk"); tfrfme.read(reinterpret_cast<char *>(&exri), sizeof(double)); exri = tfrfme->get_param("exri"); tfrfme.read(reinterpret_cast<char *>(&an), sizeof(double)); an = tfrfme->get_param("an"); tfrfme.read(reinterpret_cast<char *>(&ff), sizeof(double)); ff = tfrfme->get_param("ff"); tfrfme.read(reinterpret_cast<char *>(&tra), sizeof(double)); tra = tfrfme->get_param("tra"); tfrfme.read(reinterpret_cast<char *>(&spd), sizeof(double)); spd = tfrfme->get_param("spd"); tfrfme.read(reinterpret_cast<char *>(&frsh), sizeof(double)); frsh = tfrfme->get_param("frsh"); tfrfme.read(reinterpret_cast<char *>(&exril), sizeof(double)); exril = tfrfme->get_param("exril"); vkv = new double[nkv](); vkv = new double[nkv](); xv = new double[nxv](); yv = new double[nyv](); zv = new double[nzv](); for (int xi = 0; xi < nxv; xi++) tfrfme.read(reinterpret_cast<char *>(&(xv[xi])), sizeof(double)); for (int yi = 0; yi < nxv; yi++) tfrfme.read(reinterpret_cast<char *>(&(yv[yi])), sizeof(double)); for (int zi = 0; zi < nxv; zi++) tfrfme.read(reinterpret_cast<char *>(&(zv[zi])), sizeof(double)); fstream temptape2; fstream temptape2; string tempname2 = output_path + "/c_TEMPTAPE2"; string tempname2 = output_path + "/c_TEMPTAPE2"; temptape2.open(tempname2.c_str(), ios::in | ios::binary); temptape2.open(tempname2.c_str(), ios::in | ios::binary); Loading Loading @@ -116,15 +119,6 @@ void frfme(string data_file, string output_path) { } else { } else { printf("ERROR: could not open TEMPTAPE2 file.\n"); printf("ERROR: could not open TEMPTAPE2 file.\n"); } } for (int ixyz12 = 0; ixyz12 < nrvc; ixyz12++) { for (int j12 = 0; j12 < jlmf - 1; j12++) { double vreal, vimag; tfrfme.read(reinterpret_cast<char *>(&vreal), sizeof(double)); tfrfme.read(reinterpret_cast<char *>(&vimag), sizeof(double)); wsum[j12][ixyz12] = complex<double>(vreal, vimag); } // j12 loop } // ixyz12 loop tfrfme.close(); } else { } else { printf("ERROR: could not open TFRFME file.\n"); printf("ERROR: could not open TFRFME file.\n"); } } Loading Loading @@ -179,23 +173,22 @@ void frfme(string data_file, string output_path) { re = regex("[0-9]+"); re = regex("[0-9]+"); regex_search(str_target, m, re); regex_search(str_target, m, re); int ixi = stoi(m.str()); int ixi = stoi(m.str()); fstream tedf; string tedf_name = output_path + "/" + namef + ".hd5"; string tedf_name = output_path + "/" + namef; ScattererConfiguration *tedf = ScattererConfiguration::from_binary(tedf_name, "HDF5"); tedf.open(tedf_name.c_str(), ios::in | ios::binary); if (tedf != NULL) { if (tedf.is_open()) { int iduml, idum; int iduml, idum; tedf.read(reinterpret_cast<char *>(&iduml), sizeof(int)); iduml = (int)tedf->get_param("nsph"); for (int i = 0; i < iduml; i++) tedf.read(reinterpret_cast<char *>(&idum), sizeof(int)); idum = tedf->get_iog(iduml - 1); tedf.read(reinterpret_cast<char *>(&exdc), sizeof(double)); exdc = tedf->get_param("exdc"); tedf.read(reinterpret_cast<char *>(&wp), sizeof(double)); wp = tedf->get_param("wp"); tedf.read(reinterpret_cast<char *>(&xip), sizeof(double)); xip = tedf->get_param("xip"); tedf.read(reinterpret_cast<char *>(&idfc), sizeof(int)); idfc = (int)tedf->get_param("idfc"); tedf.read(reinterpret_cast<char *>(&nxi), sizeof(int)); nxi = (int)tedf->get_param("nxi"); if (idfc >= 0) { if (idfc >= 0) { if (ixi <= nxi) { if (ixi <= nxi) { for (int i = 0; i < ixi; i++) tedf.read(reinterpret_cast<char *>(&xi), sizeof(double)); xi = tedf->get_scale(ixi - 1); } else { // label 96 } else { // label 96 tedf.close(); delete tedf; // label 98 // label 98 string output_name = output_path + "/c_OFRFME"; string output_name = output_path + "/c_OFRFME"; FILE *output = fopen(output_name.c_str(), "w"); FILE *output = fopen(output_name.c_str(), "w"); Loading @@ -206,7 +199,7 @@ void frfme(string data_file, string output_path) { xi = xip; xi = xip; } } // label 20 // label 20 tedf.close(); delete tedf; double wn = wp / 3.0e8; double wn = wp / 3.0e8; vk = xi * wn; vk = xi * wn; exri = sqrt(exdc); exri = sqrt(exdc); Loading Loading @@ -244,26 +237,24 @@ void frfme(string data_file, string output_path) { int nxs = nxsh * 2; int nxs = nxsh * 2; int nxv = nxs + 1; int nxv = nxs + 1; int nxshpo = nxsh + 1; int nxshpo = nxsh + 1; xv = new double[nxv](); for (int i24 = nxshpo; i24 <= nxs; i24++) { xv[i24] = xv[i24 - 1] + delxyz; xv[nxv - i24 - 1] = -xv[i24]; } // i24 loop int nys = nysh * 2; int nys = nysh * 2; int nyv = nys + 1; int nyv = nys + 1; int nyshpo = nysh + 1; int nyshpo = nysh + 1; yv = new double[nyv](); for (int i25 = nyshpo; i25 <= nys; i25++) { yv[i25] = yv[i25 - 1] + delxyz; yv[nyv - i25 - 1] = -yv[i25]; } // i25 loop int nzs = nzsh * 2; int nzs = nzsh * 2; int nzv = nzs + 1; int nzv = nzs + 1; int nzshpo = nzsh + 1; int nzshpo = nzsh + 1; zv = new double[nzv](); tfrfme = new TFRFME(lmode, lm, nkv, nxv, nyv, nzv); for (int i24 = nxshpo; i24 <= nxs; i24++) { tfrfme->set_x(i24, tfrfme->get_x(i24 - 1) + delxyz); tfrfme->set_x(nxv - i24 - 1, -tfrfme->get_x(i24)); } // i24 loop for (int i25 = nyshpo; i25 <= nys; i25++) { tfrfme->set_y(i25, tfrfme->get_y(i25 - 1) + delxyz); tfrfme->set_y(nyv - i25 - 1, -tfrfme->get_y(i25)); } // i25 loop for (int i27 = nzshpo; i27 <= nzs; i27++) { for (int i27 = nzshpo; i27 <= nzs; i27++) { zv[i27] = zv[i27 - 1] + delxyz; tfrfme->set_z(i27, tfrfme->get_z(i27 - 1) + delxyz); zv[nzv - i27 - 1] = -zv[i27]; tfrfme->set_z(nzv - i27 - 1, -tfrfme->get_z(i27)); } // i27 loop } // i27 loop int nrvc = nxv * nyv * nzv; int nrvc = nxv * nyv * nzv; int nkshpo = nksh + 1; int nkshpo = nksh + 1; Loading @@ -271,28 +262,15 @@ void frfme(string data_file, string output_path) { vkv[i28] = vkv[i28 - 1] + delk; vkv[i28] = vkv[i28 - 1] + delk; vkv[nkv - i28 - 1] = -vkv[i28]; vkv[nkv - i28 - 1] = -vkv[i28]; } // i28 loop } // i28 loop tfrfme.open(tfrfme_name.c_str(), ios::out | ios::binary); if (tfrfme != NULL) { if (tfrfme.is_open()) { tfrfme->set_param("vk", vk); tfrfme.write(reinterpret_cast<char *>(&lmode), sizeof(int)); tfrfme->set_param("exri", exri); tfrfme.write(reinterpret_cast<char *>(&lm), sizeof(int)); tfrfme->set_param("an", an); tfrfme.write(reinterpret_cast<char *>(&nkv), sizeof(int)); tfrfme->set_param("ff", ff); tfrfme.write(reinterpret_cast<char *>(&nxv), sizeof(int)); tfrfme->set_param("tra", tra); tfrfme.write(reinterpret_cast<char *>(&nyv), sizeof(int)); tfrfme->set_param("spd", spd); tfrfme.write(reinterpret_cast<char *>(&nzv), sizeof(int)); tfrfme->set_param("frsh", frsh); tfrfme.write(reinterpret_cast<char *>(&vk), sizeof(double)); tfrfme->set_param("exril", exril); tfrfme.write(reinterpret_cast<char *>(&exri), sizeof(double)); tfrfme.write(reinterpret_cast<char *>(&an), sizeof(double)); tfrfme.write(reinterpret_cast<char *>(&ff), sizeof(double)); tfrfme.write(reinterpret_cast<char *>(&tra), sizeof(double)); tfrfme.write(reinterpret_cast<char *>(&spd), sizeof(double)); tfrfme.write(reinterpret_cast<char *>(&frsh), sizeof(double)); tfrfme.write(reinterpret_cast<char *>(&exril), sizeof(double)); for (int xi = 0; xi < nxv; xi++) tfrfme.write(reinterpret_cast<char *>(&(xv[xi])), sizeof(double)); for (int yi = 0; yi < nyv; yi++) tfrfme.write(reinterpret_cast<char *>(&(yv[yi])), sizeof(double)); for (int zi = 0; zi < nzv; zi++) tfrfme.write(reinterpret_cast<char *>(&(zv[zi])), sizeof(double)); fstream temptape1, temptape2; fstream temptape1, temptape2; string temp_name1 = output_path + "/c_TEMPTAPE1"; string temp_name1 = output_path + "/c_TEMPTAPE1"; string temp_name2 = output_path + "/c_TEMPTAPE2"; string temp_name2 = output_path + "/c_TEMPTAPE2"; Loading Loading @@ -330,11 +308,6 @@ void frfme(string data_file, string output_path) { } } } // jy40 loop } // jy40 loop temptape2.close(); temptape2.close(); if (wsum != NULL) { for (int wsi = wsum_size - 1; wsi > -1; wsi--) delete[] wsum[wsi]; delete[] wsum; } wsum = new complex<double>*[nlmmt]; for (int j80 = jlmf - 1; j80 < jlml; j80++) { for (int j80 = jlmf - 1; j80 < jlml; j80++) { // w matrix // w matrix if (w != NULL) { if (w != NULL) { Loading @@ -359,13 +332,13 @@ void frfme(string data_file, string output_path) { } // jy50 loop } // jy50 loop temptape1.close(); temptape1.close(); int ixyz = 0; int ixyz = 0; wsum[j80] = new complex<double>[nrvc](); for (int wj = 0; wj < nrvc; wj++) tfrfme->set_matrix_element(j80, wj, cc0); for (int iz75 = 0; iz75 < nzv; iz75++) { for (int iz75 = 0; iz75 < nzv; iz75++) { double z = zv[iz75] + frsh; double z = tfrfme->get_z(iz75) + frsh; for (int iy70 = 0; iy70 < nyv; iy70++) { for (int iy70 = 0; iy70 < nyv; iy70++) { double y = yv[iy70]; double y = tfrfme->get_y(iy70); for (int ix65 = 0; ix65 < nxv; ix65++) { for (int ix65 = 0; ix65 < nxv; ix65++) { double x = xv[ix65]; double x = tfrfme->get_x(ix65); ixyz++; ixyz++; complex<double> sumy = cc0; complex<double> sumy = cc0; for (int jy60 = 0; jy60 < nkv; jy60++) { for (int jy60 = 0; jy60 < nkv; jy60++) { Loading @@ -385,41 +358,29 @@ void frfme(string data_file, string output_path) { if (jy60 == 0 || jy60 == nkv - 1) sumx *= 0.5; if (jy60 == 0 || jy60 == nkv - 1) sumx *= 0.5; sumy += sumx; sumy += sumx; } // jy60 loop } // jy60 loop wsum[j80][ixyz - 1] = sumy * delks; tfrfme->set_matrix_element(j80, ixyz - 1, sumy * delks); } // ix65 loop } // ix65 loop } // iy70 loop } // iy70 loop } // iz75 loop } // iz75 loop } // j80 loop } // j80 loop if (jlmf != 1) { if (jlmf != 1) { tfrfme.open(tfrfme_name, ios::in | ios::out | ios::binary); lmode = (int)tfrfme->get_param("lmode"); tfrfme.read(reinterpret_cast<char *>(&lmode), sizeof(int)); lm = (int)tfrfme->get_param("lm"); tfrfme.read(reinterpret_cast<char *>(&lm), sizeof(int)); nkv = (int)tfrfme->get_param("nkv"); tfrfme.read(reinterpret_cast<char *>(&nkv), sizeof(int)); nxv = (int)tfrfme->get_param("nxv"); tfrfme.read(reinterpret_cast<char *>(&nxv), sizeof(int)); nyv = (int)tfrfme->get_param("nyv"); tfrfme.read(reinterpret_cast<char *>(&nyv), sizeof(int)); nzv = (int)tfrfme->get_param("nzv"); tfrfme.read(reinterpret_cast<char *>(&nzv), sizeof(int)); vk = tfrfme->get_param("vk"); tfrfme.read(reinterpret_cast<char *>(&vk), sizeof(double)); exri = tfrfme->get_param("exri"); tfrfme.read(reinterpret_cast<char *>(&exri), sizeof(double)); an = tfrfme->get_param("an"); tfrfme.read(reinterpret_cast<char *>(&an), sizeof(double)); ff = tfrfme->get_param("ff"); tfrfme.read(reinterpret_cast<char *>(&ff), sizeof(double)); tra = tfrfme->get_param("tra"); tfrfme.read(reinterpret_cast<char *>(&tra), sizeof(double)); spd = tfrfme->get_param("spd"); tfrfme.read(reinterpret_cast<char *>(&spd), sizeof(double)); frsh = tfrfme->get_param("frsh"); tfrfme.read(reinterpret_cast<char *>(&frsh), sizeof(double)); exril = tfrfme->get_param("exril"); tfrfme.read(reinterpret_cast<char *>(&exril), sizeof(double)); for (int i = 0; i < nxv; i++) tfrfme.read(reinterpret_cast<char *>(&(xv[i])), sizeof(double)); for (int i = 0; i < nyv; i++) tfrfme.read(reinterpret_cast<char *>(&(yv[i])), sizeof(double)); for (int i = 0; i < nzv; i++) tfrfme.read(reinterpret_cast<char *>(&(zv[i])), sizeof(double)); } } // label 88 // label 88 for (int ixyz = 0; ixyz < nrvc; ixyz++) { tfrfme->write_binary(tfrfme_name, "HDF5"); for (int j = 0; j < jlml; j++) { double vreal = wsum[j][ixyz].real(); double vimag = wsum[j][ixyz].imag(); tfrfme.write(reinterpret_cast<char *>(&vreal), sizeof(double)); tfrfme.write(reinterpret_cast<char *>(&vimag), sizeof(double)); } // j loop } // ixyz loop tfrfme.close(); string output_name = output_path + "/c_OFRFME"; string output_name = output_path + "/c_OFRFME"; FILE *output = fopen(output_name.c_str(), "w"); FILE *output = fopen(output_name.c_str(), "w"); fprintf(output, " IF JLML < NLMMT, PRESERVE TEMPTAPE1, TEMPTAPE2, AND TFRFRME,\n"); fprintf(output, " IF JLML < NLMMT, PRESERVE TEMPTAPE1, TEMPTAPE2, AND TFRFRME,\n"); Loading @@ -441,11 +402,8 @@ void frfme(string data_file, string output_path) { } } } } // label 45 // label 45 if (tfrfme.is_open()) tfrfme.close(); if (tfrfme != NULL) delete tfrfme; delete[] file_lines; delete[] file_lines; if (xv != NULL) delete[] xv; if (yv != NULL) delete[] yv; if (zv != NULL) delete[] zv; if (vkv != NULL) delete[] vkv; if (vkv != NULL) delete[] vkv; if (vkzm != NULL) { if (vkzm != NULL) { for (int vki = nkv - 1; vki > -1; vki--) delete[] vkzm[vki]; for (int vki = nkv - 1; vki > -1; vki--) delete[] vkzm[vki]; Loading @@ -455,10 +413,6 @@ void frfme(string data_file, string output_path) { for (int wi = nkv - 1; wi > -1; wi--) delete[] w[wi]; for (int wi = nkv - 1; wi > -1; wi--) delete[] w[wi]; delete[] w; delete[] w; } } if (wsum != NULL) { for (int wsi = wsum_size - 1; wsi > -1; wsi--) delete[] wsum[wsi]; delete[] wsum; } if (wk != NULL) delete[] wk; if (wk != NULL) delete[] wk; printf("Done.\n"); printf("Done.\n"); } }
src/trapping/clffft.cpp +58 −46 Original line number Original line Diff line number Diff line Loading @@ -4,6 +4,7 @@ */ */ #include <complex> #include <complex> #include <cstdio> #include <cstdio> #include <exception> #include <fstream> #include <fstream> #include <regex> #include <regex> #include <string> #include <string> Loading @@ -20,6 +21,10 @@ #include "../include/sph_subs.h" #include "../include/sph_subs.h" #endif #endif #ifndef INCLUDE_TFRFME_H_ #include "../include/tfrfme.h" #endif #ifndef INCLUDE_TRA_SUBS_H_ #ifndef INCLUDE_TRA_SUBS_H_ #include "../include/tra_subs.h" #include "../include/tra_subs.h" #endif #endif Loading @@ -38,7 +43,6 @@ void lffft(string data_file, string output_path) { fstream tlfff, tlfft; fstream tlfff, tlfft; double ****zpv = NULL; double ****zpv = NULL; double *xv = NULL, *yv = NULL, *zv = NULL; complex<double> *ac = NULL, *ws = NULL, *wsl = NULL; complex<double> *ac = NULL, *ws = NULL, *wsl = NULL; complex<double> **am0m = NULL; complex<double> **am0m = NULL; complex<double> **amd = NULL; complex<double> **amd = NULL; Loading Loading @@ -140,37 +144,31 @@ void lffft(string data_file, string output_path) { } } // label 150 // label 150 ttms.close(); ttms.close(); fstream binary_input; TFRFME *tfrfme = NULL; string binary_name; string binary_name; if (jss != 1) binary_name = output_path + "/c_TFRFME"; if (jss != 1) binary_name = output_path + "/c_TFRFME.hd5"; else binary_name = output_path + "/c_TWS"; else binary_name = output_path + "/c_TWS"; binary_input.open(binary_name, ios::in | ios::binary); tfrfme = TFRFME::from_binary(binary_name, "HDF5"); if (binary_input.is_open()) { if (tfrfme != NULL) { int lmode, lm, nkv, nxv, nyv, nzv; int lmode, lm, nkv, nxv, nyv, nzv; double vk, exri, an, ff, tra; double vk, exri, an, ff, tra; double spd, frsh, exril; double spd, frsh, exril; binary_input.read(reinterpret_cast<char *>(&lmode), sizeof(int)); lmode = (int)tfrfme->get_param("lmode"); binary_input.read(reinterpret_cast<char *>(&lm), sizeof(int)); lm = (int)tfrfme->get_param("lm"); binary_input.read(reinterpret_cast<char *>(&nkv), sizeof(int)); nkv = (int)tfrfme->get_param("nkv"); binary_input.read(reinterpret_cast<char *>(&nxv), sizeof(int)); nxv = (int)tfrfme->get_param("nxv"); binary_input.read(reinterpret_cast<char *>(&nyv), sizeof(int)); nyv = (int)tfrfme->get_param("nyv"); binary_input.read(reinterpret_cast<char *>(&nzv), sizeof(int)); nzv = (int)tfrfme->get_param("nzv"); if (lm >= le) { if (lm >= le) { binary_input.read(reinterpret_cast<char *>(&vk), sizeof(double)); vk = tfrfme->get_param("vk"); binary_input.read(reinterpret_cast<char *>(&exri), sizeof(double)); exri = tfrfme->get_param("exri"); binary_input.read(reinterpret_cast<char *>(&an), sizeof(double)); an = tfrfme->get_param("an"); binary_input.read(reinterpret_cast<char *>(&ff), sizeof(double)); ff = tfrfme->get_param("ff"); binary_input.read(reinterpret_cast<char *>(&tra), sizeof(double)); tra = tfrfme->get_param("tra"); if (vk == vks && exri == exris) { if (vk == vks && exri == exris) { binary_input.read(reinterpret_cast<char *>(&spd), sizeof(double)); spd = tfrfme->get_param("spd"); binary_input.read(reinterpret_cast<char *>(&frsh), sizeof(double)); frsh = tfrfme->get_param("frsh"); binary_input.read(reinterpret_cast<char *>(&exril), sizeof(double)); exril = tfrfme->get_param("exril"); xv = new double[nxv]; for (int i = 0; i < nxv; i++) binary_input.read(reinterpret_cast<char *>(&(xv[i])), sizeof(double)); yv = new double[nyv]; for (int i = 0; i < nyv; i++) binary_input.read(reinterpret_cast<char *>(&(yv[i])), sizeof(double)); zv = new double[nzv]; for (int i = 0; i < nzv; i++) binary_input.read(reinterpret_cast<char *>(&(zv[i])), sizeof(double)); bool goto160 = false; bool goto160 = false; if (jft <= 0) { if (jft <= 0) { zpv = new double***[le]; zpv = new double***[le]; Loading Loading @@ -204,12 +202,18 @@ void lffft(string data_file, string output_path) { tlfff.write(reinterpret_cast<char *>(&spd), sizeof(double)); tlfff.write(reinterpret_cast<char *>(&spd), sizeof(double)); tlfff.write(reinterpret_cast<char *>(&frsh), sizeof(double)); tlfff.write(reinterpret_cast<char *>(&frsh), sizeof(double)); tlfff.write(reinterpret_cast<char *>(&exril), sizeof(double)); tlfff.write(reinterpret_cast<char *>(&exril), sizeof(double)); for (int i = 0; i < nxv; i++) for (int i = 0; i < nxv; i++) { tlfff.write(reinterpret_cast<char *>(&(xv[i])), sizeof(double)); double x = tfrfme->get_x(i); for (int i = 0; i < nyv; i++) tlfff.write(reinterpret_cast<char *>(&x), sizeof(double)); tlfff.write(reinterpret_cast<char *>(&(yv[i])), sizeof(double)); } for (int i = 0; i < nzv; i++) for (int i = 0; i < nyv; i++) { tlfff.write(reinterpret_cast<char *>(&(zv[i])), sizeof(double)); double y = tfrfme->get_y(i); tlfff.write(reinterpret_cast<char *>(&y), sizeof(double)); } for (int i = 0; i < nzv; i++) { double z = tfrfme->get_z(i); tlfff.write(reinterpret_cast<char *>(&z), sizeof(double)); } if (jft < 0) goto160 = true; if (jft < 0) goto160 = true; } else { // Should never happen. } else { // Should never happen. printf("ERROR: could not open TLFFF file.\n"); printf("ERROR: could not open TLFFF file.\n"); Loading @@ -236,12 +240,18 @@ void lffft(string data_file, string output_path) { tlfft.write(reinterpret_cast<char *>(&spd), sizeof(double)); tlfft.write(reinterpret_cast<char *>(&spd), sizeof(double)); tlfft.write(reinterpret_cast<char *>(&frsh), sizeof(double)); tlfft.write(reinterpret_cast<char *>(&frsh), sizeof(double)); tlfft.write(reinterpret_cast<char *>(&exril), sizeof(double)); tlfft.write(reinterpret_cast<char *>(&exril), sizeof(double)); for (int i = 0; i < nxv; i++) for (int i = 0; i < nxv; i++) { tlfft.write(reinterpret_cast<char *>(&(xv[i])), sizeof(double)); double x = tfrfme->get_x(i); for (int i = 0; i < nyv; i++) tlfft.write(reinterpret_cast<char *>(&x), sizeof(double)); tlfft.write(reinterpret_cast<char *>(&(yv[i])), sizeof(double)); } for (int i = 0; i < nzv; i++) for (int i = 0; i < nyv; i++) { tlfft.write(reinterpret_cast<char *>(&(zv[i])), sizeof(double)); double y = tfrfme->get_y(i); tlfft.write(reinterpret_cast<char *>(&y), sizeof(double)); } for (int i = 0; i < nzv; i++) { double z = tfrfme->get_z(i); tlfft.write(reinterpret_cast<char *>(&z), sizeof(double)); } } else { // Should never happen. } else { // Should never happen. printf("ERROR: could not open TLFFT file.\n"); printf("ERROR: could not open TLFFT file.\n"); } } Loading @@ -262,13 +272,18 @@ void lffft(string data_file, string output_path) { for (int iy475 = 0; iy475 < nyv; iy475++) { for (int iy475 = 0; iy475 < nyv; iy475++) { for (int ix475 = 0; ix475 < nxv; ix475++) { for (int ix475 = 0; ix475 < nxv; ix475++) { for (int i = 0; i < nlmmt; i++) { for (int i = 0; i < nlmmt; i++) { double vreal, vimag; //double vreal, vimag; binary_input.read(reinterpret_cast<char *>(&vreal), sizeof(double)); //binary_input.read(reinterpret_cast<char *>(&vreal), sizeof(double)); binary_input.read(reinterpret_cast<char *>(&vimag), sizeof(double)); //binary_input.read(reinterpret_cast<char *>(&vimag), sizeof(double)); int row = i; int col = (nzv * iz475) + (nyv * iy475) + ix475; complex<double> value = tfrfme->get_matrix_element(row, col); if (lm <= le) { if (lm <= le) { ws[i] = complex<double>(vreal, vimag); //ws[i] = complex<double>(vreal, vimag); ws[i] = value; } else { // label 170 } else { // label 170 wsl[i] = complex<double>(vreal, vimag); //wsl[i] = complex<double>(vreal, vimag); wsl[i] = value; for (int i175 = 0; i175 < cil->nlem; i175++) { for (int i175 = 0; i175 < cil->nlem; i175++) { int ie = i175 + cil->nlem; int ie = i175 + cil->nlem; int iel = i175 + nlmm; int iel = i175 + nlmm; Loading Loading @@ -384,7 +399,7 @@ void lffft(string data_file, string output_path) { fclose(output67); fclose(output67); } } } } binary_input.close(); delete tfrfme; } else { } else { printf("ERROR: could not open binary input file %s.\n", binary_name.c_str()); printf("ERROR: could not open binary input file %s.\n", binary_name.c_str()); } } Loading @@ -395,9 +410,6 @@ void lffft(string data_file, string output_path) { // Clean up memory // Clean up memory if (ac != NULL) delete[] ac; if (ac != NULL) delete[] ac; if (ws != NULL) delete[] ws; if (ws != NULL) delete[] ws; if (xv != NULL) delete[] xv; if (yv != NULL) delete[] yv; if (zv != NULL) delete[] zv; if (wsl != NULL) delete[] wsl; if (wsl != NULL) delete[] wsl; if (tmsm != NULL) delete[] tmsm; if (tmsm != NULL) delete[] tmsm; if (tmse != NULL) delete[] tmse; if (tmse != NULL) delete[] tmse; Loading