diff --git a/src/cluster/cluster.cpp b/src/cluster/cluster.cpp index 06410a2b71fd03a834806a6421f2bfb39ad97025..44aff3db320a1151c8988aba22a25a5fdff428c9 100644 --- a/src/cluster/cluster.cpp +++ b/src/cluster/cluster.cpp @@ -4,14 +4,14 @@ you may not use this file except in compliance with the License. You may obtain a copy of the License at - http://www.apache.org/licenses/LICENSE-2.0 + http://www.apache.org/licenses/LICENSE-2.0 Unless required by applicable law or agreed to in writing, software distributed under the License is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. See the License for the specific language governing permissions and limitations under the License. - */ +*/ /*! \file cluster.cpp * @@ -21,6 +21,7 @@ #include <cstdio> #include <exception> #include <fstream> +#include <hdf5.h> #include <string> #ifdef _OPENMP #include <omp.h> @@ -73,11 +74,19 @@ #include "../include/algebraic.h" #endif +#ifndef INCLUDE_LIST_H_ +#include "../include/List.h" +#endif + +#ifndef INCLUDE_FILE_IO_H_ +#include "../include/file_io.h" +#endif + using namespace std; // I would like to put it all in a struct, but then I'd have to write a constructor for it, due to members defined as references, creating a worse nightmare than the one I'd like to simplify... -int cluster_jxi488_cycle(int jxi488, ScattererConfiguration *sconf, GeometryConfiguration *gconf, ScatteringAngles *sa, ClusterIterationData *cid, FILE *output, const string& output_path, fstream& tppoan); +int cluster_jxi488_cycle(int jxi488, ScattererConfiguration *sconf, GeometryConfiguration *gconf, ScatteringAngles *sa, ClusterIterationData *cid, VirtualAsciiFile *output, const string& output_path, VirtualBinaryFile *vtppoanp); /*! \brief C++ implementation of CLU * @@ -152,54 +161,80 @@ void cluster(const string& config_file, const string& data_file, const string& o // Shortcuts to variables stored in configuration objects ScatteringAngles *p_scattering_angles = new ScatteringAngles(gconf); double wp = sconf->wp; - FILE *output = fopen((output_path + "/c_OCLU").c_str(), "w"); + //FILE *output = fopen((output_path + "/c_OCLU").c_str(), "w"); + VirtualAsciiFile *p_output = new VirtualAsciiFile(); + // for the time being, this is ok. When we can, add some logic in the sprintf calls that checks if a longer buffer would be needed, and in case expand it + // in any case, replace all sprintf() with snprintf(), to avoid in any case writing more than the available buffer size + char virtual_line[256]; ClusterIterationData *cid = new ClusterIterationData(gconf, sconf, mpidata); const int ndi = cid->c4->nsph * cid->c4->nlim; np_int ndit = 2 * ndi; logger->log("INFO: Size of matrices to invert: " + to_string((int64_t)ndit) + " x " + to_string((int64_t)ndit) +".\n"); time_logger->log("INFO: Size of matrices to invert: " + to_string((int64_t)ndit) + " x " + to_string((int64_t)ndit) +".\n"); - fprintf(output, " READ(IR,*)NSPH,LI,LE,MXNDM,INPOL,NPNT,NPNTTS,IAVM,ISAM\n"); - fprintf(output, " %5d%5d%5d%5ld%5d%5d%5d%5d%5d\n", + sprintf(virtual_line, " READ(IR,*)NSPH,LI,LE,MXNDM,INPOL,NPNT,NPNTTS,IAVM,ISAM\n\0"); + p_output->append_line(virtual_line); + sprintf(virtual_line, " %5d%5d%5d%5ld%5d%5d%5d%5d%5d\n\0", nsph, cid->c4->li, cid->c4->le, gconf->mxndm, gconf->in_pol, gconf->npnt, gconf->npntts, gconf->iavm, gconf->iavm ); - fprintf(output, " READ(IR,*)RXX(I),RYY(I),RZZ(I)\n"); - for (int ri = 0; ri < nsph; ri++) fprintf(output, "%17.8lE%17.8lE%17.8lE\n", - gconf->get_sph_x(ri), gconf->get_sph_y(ri), gconf->get_sph_z(ri) - ); - fprintf(output, " READ(IR,*)TH,THSTP,THLST,THS,THSSTP,THSLST\n"); - fprintf( - output, " %10.3lE%10.3lE%10.3lE%10.3lE%10.3lE%10.3lE\n", + p_output->append_line(virtual_line); + sprintf(virtual_line, " READ(IR,*)RXX(I),RYY(I),RZZ(I)\n\0"); + p_output->append_line(virtual_line); + for (int ri = 0; ri < nsph; ri++) { + sprintf(virtual_line, "%17.8lE%17.8lE%17.8lE\n\0", + gconf->get_sph_x(ri), gconf->get_sph_y(ri), gconf->get_sph_z(ri) + ); + p_output->append_line(virtual_line); + } + sprintf(virtual_line, " READ(IR,*)TH,THSTP,THLST,THS,THSSTP,THSLST\n\0"); + p_output->append_line(virtual_line); + sprintf( + virtual_line, " %10.3lE%10.3lE%10.3lE%10.3lE%10.3lE%10.3lE\n\0", p_scattering_angles->th, p_scattering_angles->thstp, p_scattering_angles->thlst, p_scattering_angles->ths, p_scattering_angles->thsstp, p_scattering_angles->thslst ); - fprintf(output, " READ(IR,*)PH,PHSTP,PHLST,PHS,PHSSTP,PHSLST\n"); - fprintf( - output, " %10.3lE%10.3lE%10.3lE%10.3lE%10.3lE%10.3lE\n", + p_output->append_line(virtual_line); + sprintf(virtual_line, " READ(IR,*)PH,PHSTP,PHLST,PHS,PHSSTP,PHSLST\n\0"); + p_output->append_line(virtual_line); + sprintf( + virtual_line, " %10.3lE%10.3lE%10.3lE%10.3lE%10.3lE%10.3lE\n\0", p_scattering_angles->ph, p_scattering_angles->phstp, p_scattering_angles->phlst, p_scattering_angles->phs, p_scattering_angles->phsstp, p_scattering_angles->phslst ); - fprintf(output, " READ(IR,*)JWTM\n"); - fprintf(output, " %5d\n", gconf->jwtm); - fprintf(output, " READ(ITIN)NSPHT\n"); - fprintf(output, " READ(ITIN)(IOG(I),I=1,NSPH)\n"); - fprintf(output, " READ(ITIN)EXDC,WP,XIP,IDFC,NXI\n"); - fprintf(output, " READ(ITIN)(XIV(I),I=1,NXI)\n"); - fprintf(output, " READ(ITIN)NSHL(I),ROS(I)\n"); - fprintf(output, " READ(ITIN)(RCF(I,NS),NS=1,NSH)\n"); - fprintf(output, " \n"); + p_output->append_line(virtual_line); + sprintf(virtual_line, " READ(IR,*)JWTM\n\0"); + p_output->append_line(virtual_line); + sprintf(virtual_line, " %5d\n\0", gconf->jwtm); + p_output->append_line(virtual_line); + sprintf(virtual_line, " READ(ITIN)NSPHT\n\0"); + p_output->append_line(virtual_line); + sprintf(virtual_line, " READ(ITIN)(IOG(I),I=1,NSPH)\n\0"); + p_output->append_line(virtual_line); + sprintf(virtual_line, " READ(ITIN)EXDC,WP,XIP,IDFC,NXI\n\0"); + p_output->append_line(virtual_line); + sprintf(virtual_line, " READ(ITIN)(XIV(I),I=1,NXI)\n\0"); + p_output->append_line(virtual_line); + sprintf(virtual_line, " READ(ITIN)NSHL(I),ROS(I)\n\0"); + p_output->append_line(virtual_line); + sprintf(virtual_line, " READ(ITIN)(RCF(I,NS),NS=1,NSH)\n\0"); + p_output->append_line(virtual_line); + sprintf(virtual_line, " \n\0"); + p_output->append_line(virtual_line); str(sconf, cid->c1, cid->c1ao, cid->c3, cid->c4, cid->c6); thdps(cid->c4->lm, cid->zpv); double exdc = sconf->exdc; double exri = sqrt(exdc); - fprintf(output, " REFR. INDEX OF EXTERNAL MEDIUM=%15.7lE\n", exri); - fstream *tppoanp = new fstream; - fstream &tppoan = *tppoanp; + sprintf(virtual_line, " REFR. INDEX OF EXTERNAL MEDIUM=%15.7lE\n\0", exri); + p_output->append_line(virtual_line); + VirtualBinaryFile *vtppoanp = new VirtualBinaryFile(); + VirtualBinaryLine *vbinline = NULL; + // fstream *tppoanp = new fstream; + // fstream &tppoan = *tppoanp; string tppoan_name = output_path + "/c_TPPOAN"; - tppoan.open(tppoan_name.c_str(), ios::out | ios::binary); - if (tppoan.is_open()) { + // tppoan.open(tppoan_name.c_str(), ios::out | ios::binary); + // if (tppoan.is_open()) { #ifdef USE_MAGMA logger->log("INFO: using MAGMA calls.\n", LOG_INFO); #elif defined USE_LAPACK @@ -215,18 +250,28 @@ void cluster(const string& config_file, const string& data_file, const string& o int nths = p_scattering_angles->nths; int nph = p_scattering_angles->nph; int nphs = p_scattering_angles->nphs; - tppoan.write(reinterpret_cast<char *>(&iavm), sizeof(int)); - tppoan.write(reinterpret_cast<char *>(&isam), sizeof(int)); - tppoan.write(reinterpret_cast<char *>(&inpol), sizeof(int)); - tppoan.write(reinterpret_cast<char *>(&nxi), sizeof(int)); - tppoan.write(reinterpret_cast<char *>(&nth), sizeof(int)); - tppoan.write(reinterpret_cast<char *>(&nph), sizeof(int)); - tppoan.write(reinterpret_cast<char *>(&nths), sizeof(int)); - tppoan.write(reinterpret_cast<char *>(&nphs), sizeof(int)); + // tppoan.write(reinterpret_cast<char *>(&iavm), sizeof(int)); + vtppoanp->append_line(VirtualBinaryLine(iavm)); + // tppoan.write(reinterpret_cast<char *>(&isam), sizeof(int)); + vtppoanp->append_line(VirtualBinaryLine(isam)); + // tppoan.write(reinterpret_cast<char *>(&inpol), sizeof(int)); + vtppoanp->append_line(VirtualBinaryLine(inpol)); + // tppoan.write(reinterpret_cast<char *>(&nxi), sizeof(int)); + vtppoanp->append_line(VirtualBinaryLine(nxi)); + // tppoan.write(reinterpret_cast<char *>(&nth), sizeof(int)); + vtppoanp->append_line(VirtualBinaryLine(nth)); + // tppoan.write(reinterpret_cast<char *>(&nph), sizeof(int)); + vtppoanp->append_line(VirtualBinaryLine(nph)); + // tppoan.write(reinterpret_cast<char *>(&nths), sizeof(int)); + vtppoanp->append_line(VirtualBinaryLine(nths)); + // tppoan.write(reinterpret_cast<char *>(&nphs), sizeof(int)); + vtppoanp->append_line(VirtualBinaryLine(nphs)); if (sconf->idfc < 0) { cid->vk = cid->xip * cid->wn; - fprintf(output, " VK=%15.7lE, XI IS SCALE FACTOR FOR LENGTHS\n", cid->vk); - fprintf(output, " \n"); + sprintf(virtual_line, " VK=%15.7lE, XI IS SCALE FACTOR FOR LENGTHS\n\0", cid->vk); + p_output->append_line(virtual_line); + sprintf(virtual_line, " \n\0"); + p_output->append_line(virtual_line); } // do the first iteration on jxi488 separately, since it seems to be different from the others int jxi488 = 1; @@ -234,7 +279,7 @@ void cluster(const string& config_file, const string& data_file, const string& o #ifdef USE_NVTX nvtxRangePush("First iteration"); #endif - int jer = cluster_jxi488_cycle(jxi488, sconf, gconf, p_scattering_angles, cid, output, output_path, tppoan); + int jer = cluster_jxi488_cycle(jxi488, sconf, gconf, p_scattering_angles, cid, p_output, output_path, vtppoanp); #ifdef USE_NVTX nvtxRangePop(); #endif @@ -249,9 +294,10 @@ void cluster(const string& config_file, const string& data_file, const string& o time_logger->log(message); if (jer != 0) { // First loop failed. Halt the calculation. - tppoan.close(); + // tppoan.close(); fclose(timing_file); - fclose(output); + //fclose(output); + delete p_output; delete p_scattering_angles; delete cid; delete logger; @@ -260,6 +306,15 @@ void cluster(const string& config_file, const string& data_file, const string& o delete gconf; return; } + // do the first outputs here, so that I open here the new files, afterwards I only append + p_output->write_to_disk(output_path + "/c_OCLU"); + // reallocate a new one (even if it would be more efficient to emty the existing one + delete p_output; + p_output = new VirtualAsciiFile(); + // now tppoan + vtppoanp->write_to_disk(output_path + "/c_TPPOAN"); + delete vtppoanp; + vtppoanp = new VirtualBinaryFile(); // here go the calls that send data to be duplicated on other MPI processes from process 0 to others, using MPI broadcasts, but only if MPI is actually used #ifdef MPI_VERSION @@ -272,6 +327,8 @@ void cluster(const string& config_file, const string& data_file, const string& o #endif // Create this variable and initialise it with a default here, so that it is defined anyway, with or without OpenMP support enabled int ompnumthreads = 1; + VirtualAsciiFile **p_outarray = NULL; + VirtualBinaryFile **vtppoanarray = NULL; #ifdef USE_NVTX nvtxRangePush("Parallel loop"); @@ -285,23 +342,35 @@ void cluster(const string& config_file, const string& data_file, const string& o myompthread = omp_get_thread_num(); if (myompthread == 0) ompnumthreads = omp_get_num_threads(); #endif + if (myompthread == 0) { + p_outarray = new VirtualAsciiFile*[ompnumthreads]; + vtppoanarray = new VirtualBinaryFile*[ompnumthreads]; + } // To test parallelism, I will now start feeding this function with "clean" copies of the parameters, so that they will not be changed by previous iterations, and each one will behave as the first one. Define all (empty) variables here, so they have the correct scope, then they get different definitions depending on thread number ClusterIterationData *cid_2 = NULL; - FILE *output_2 = NULL; - fstream *tppoanp_2 = NULL; + //FILE *output_2 = NULL; + VirtualAsciiFile *p_output_2 = NULL; + VirtualBinaryFile *vtppoanp_2 = NULL; + // fstream *tppoanp_2 = NULL; // for threads other than the 0, create distinct copies of all relevant data, while for thread 0 just define new references / pointers to the original ones if (myompthread == 0) { cid_2 = cid; - output_2 = output; - tppoanp_2 = tppoanp; + //output_2 = output; + p_output_2 = p_output; + // tppoanp_2 = tppoanp; + vtppoanp_2 = vtppoanp; } else { // this is not thread 0, so do create fresh copies of all local variables cid_2 = new ClusterIterationData(*cid); - output_2 = fopen((output_path + "/c_OCLU_" + to_string(mpidata->rank) + "_" + to_string(myompthread)).c_str(), "w"); - tppoanp_2 = new fstream; - tppoanp_2->open((output_path + "/c_TPPOAN_" + to_string(mpidata->rank) + "_" + to_string(myompthread)).c_str(), ios::out | ios::binary); + //output_2 = fopen((output_path + "/c_OCLU_" + to_string(mpidata->rank) + "_" + to_string(myompthread)).c_str(), "w"); + p_output_2 = new VirtualAsciiFile(); + vtppoanp_2 = new VirtualBinaryFile(); + // tppoanp_2 = new fstream; + // tppoanp_2->open((output_path + "/c_TPPOAN_" + to_string(mpidata->rank) + "_" + to_string(myompthread)).c_str(), ios::out | ios::binary); } - fstream &tppoan_2 = *tppoanp_2; + p_outarray[myompthread] = p_output_2; + vtppoanarray[myompthread] = vtppoanp_2; + // fstream &tppoan_2 = *tppoanp_2; // make sure all threads align here: I don't want the following loop to accidentally start for thread 0, possibly modifying some variables before they are copied by all other threads #pragma omp barrier if (myompthread==0) { @@ -310,16 +379,18 @@ void cluster(const string& config_file, const string& data_file, const string& o // ok, now I can actually start the parallel calculations #pragma omp for for (jxi488 = cid_2->firstxi; jxi488 <= cid_2->lastxi; jxi488++) { - int jer = cluster_jxi488_cycle(jxi488, sconf, gconf, p_scattering_angles, cid_2, output_2, output_path, *tppoanp_2); + int jer = cluster_jxi488_cycle(jxi488, sconf, gconf, p_scattering_angles, cid_2, p_output_2, output_path, vtppoanp_2); } #pragma omp barrier // only threads different from 0 have to free local copies of variables and close local files if (myompthread != 0) { delete cid_2; - fclose(output_2); - tppoanp_2->close(); - delete tppoanp_2; + //fclose(output_2); + // p_output_2->write_to_disk(output_path + "/c_OCLU_" + to_string(mpidata->rank) + "_" + to_string(myompthread)); + // delete p_output_2; + // tppoanp_2->close(); + // delete tppoanp_2; } #pragma omp barrier { @@ -332,88 +403,126 @@ void cluster(const string& config_file, const string& data_file, const string& o nvtxRangePush("Output concatenation"); #endif -#ifdef _OPENMP #pragma omp barrier { // thread 0 already wrote on global files, skip it and take care of appending the others - for (int ri = 1; ri < ompnumthreads; ri++) { - string partial_file_name = output_path + "/c_OCLU_" + to_string(mpidata->rank) + "_" + to_string(ri); - string message = "Copying ASCII output in MPI process " + to_string(mpidata->rank) + " of thread " + to_string(ri) + " of " + to_string(ompnumthreads - 1) + "... "; - logger->log(message, LOG_DEBG); - FILE *partial_output = fopen(partial_file_name.c_str(), "r"); - int c = fgetc(partial_output); - while (c != EOF) { - fputc(c, output); - c = fgetc(partial_output); - } - fclose(partial_output); - remove(partial_file_name.c_str()); - logger->log("done.\n", LOG_DEBG); - partial_file_name = output_path + "/c_TPPOAN_" + to_string(mpidata->rank) + "_" + to_string(ri); - message = "Copying binary output in MPI process " + to_string(mpidata->rank) + " of thread " + to_string(ri) + " of " + to_string(ompnumthreads - 1) + "... "; - logger->log(message, LOG_DEBG); - fstream partial_tppoan; - partial_tppoan.open(partial_file_name.c_str(), ios::in | ios::binary); - partial_tppoan.seekg(0, ios::end); - long buffer_size = partial_tppoan.tellg(); - char *binary_buffer = new char[buffer_size]; - partial_tppoan.seekg(0, ios::beg); - partial_tppoan.read(binary_buffer, buffer_size); - tppoan.write(binary_buffer, buffer_size); - partial_tppoan.close(); - delete[] binary_buffer; - remove(partial_file_name.c_str()); - logger->log("done.\n", LOG_DEBG); + for (int ti=1; ti<ompnumthreads; ti++) { + p_outarray[0]->append(*(p_outarray[ti])); + delete p_outarray[ti]; + vtppoanarray[0]->append(*(vtppoanarray[ti])); + delete vtppoanarray[ti]; } + p_outarray[0]->append_to_disk(output_path + "/c_OCLU"); + delete p_outarray[0]; + delete[] p_outarray; + vtppoanarray[0]->append_to_disk(output_path + "/c_TPPOAN"); + delete vtppoanarray[0]; + delete[] vtppoanarray; + // for (int ri = 1; ri < ompnumthreads; ri++) { + // string partial_file_name = output_path + "/c_OCLU_" + to_string(mpidata->rank) + "_" + to_string(ri); + // string message = "Copying ASCII output in MPI process " + to_string(mpidata->rank) + " of thread " + to_string(ri) + " of " + to_string(ompnumthreads - 1) + "... "; + // logger->log(message, LOG_DEBG); + // FILE *partial_output = fopen(partial_file_name.c_str(), "r"); + // // have a look at the getline() or fgets() C functions to read one line at a time, instead of char by char, it would simplify things + // char virtual_line[256]; + // int index = 0; + // int c = fgetc(partial_output); + // while (c != EOF) { + // virtual_line[index++] = c; + // if (c == '\n') { + // virtual_line[index] = '\0'; + // p_output->append_line(virtual_line); + // index = 0; + // } + // c = fgetc(partial_output); + // } + // fclose(partial_output); + // remove(partial_file_name.c_str()); + // logger->log("done.\n", LOG_DEBG); + // string partial_file_name = output_path + "/c_TPPOAN_" + to_string(mpidata->rank) + "_" + to_string(ri); + // string message = "Copying binary output in MPI process " + to_string(mpidata->rank) + " of thread " + to_string(ri) + " of " + to_string(ompnumthreads - 1) + "... "; + // logger->log(message, LOG_DEBG); + // fstream partial_tppoan; + // partial_tppoan.open(partial_file_name.c_str(), ios::in | ios::binary); + // partial_tppoan.seekg(0, ios::end); + // long buffer_size = partial_tppoan.tellg(); + // char *binary_buffer = new char[buffer_size]; + // partial_tppoan.seekg(0, ios::beg); + // partial_tppoan.read(binary_buffer, buffer_size); + // // tppoan.write(binary_buffer, buffer_size); + // partial_tppoan.close(); + // delete[] binary_buffer; + // remove(partial_file_name.c_str()); + // logger->log("done.\n", LOG_DEBG); + // } } -#endif // here go the code to append the files written in MPI processes > 0 to the ones on MPI process 0 #ifdef MPI_VERSION if (mpidata->mpirunning) { // only go through this if MPI has been actually used for (int rr=1; rr<mpidata->nprocs; rr++) { // get the data from process rr - // how many openmp threads did process rr use? - int remotethreads; - MPI_Recv(&remotethreads, 1, MPI_INT, rr, 1, MPI_COMM_WORLD, MPI_STATUS_IGNORE); - for (int ri=0; ri<remotethreads; ri++) { - // first get the ASCII local file - char *chunk_buffer; - int chunk_buffer_size = -1; - MPI_Recv(&chunk_buffer_size, 1, MPI_INT, rr, 1, MPI_COMM_WORLD, MPI_STATUS_IGNORE); - while (chunk_buffer_size != 0) { - char *chunk_buffer = new char[chunk_buffer_size]; - MPI_Recv(chunk_buffer, chunk_buffer_size, MPI_CHAR, rr, 0, MPI_COMM_WORLD, MPI_STATUS_IGNORE); - fputs(chunk_buffer, output); - delete[] chunk_buffer; - MPI_Recv(&chunk_buffer_size, 1, MPI_INT, rr, 1, MPI_COMM_WORLD, MPI_STATUS_IGNORE); - } - // if (ri<remotethreads-1) fprintf(output, "\n"); + VirtualAsciiFile *p_output = new VirtualAsciiFile(mpidata, rr); + p_output->append_to_disk(output_path + "/c_OCLU"); + delete p_output; + VirtualBinaryFile *vtppoanp = new VirtualBinaryFile(mpidata, rr); + vtppoanp->append_to_disk(output_path + "/c_TPPOAN_bis"); + delete vtppoanp; + // // how many openmp threads did process rr use? + // int remotethreads; + // MPI_Recv(&remotethreads, 1, MPI_INT, rr, 1, MPI_COMM_WORLD, MPI_STATUS_IGNORE); + // for (int ri=0; ri<remotethreads; ri++) { + // // first get the ASCII local file + // char *chunk_buffer; + // int chunk_buffer_size = -1; + // MPI_Recv(&chunk_buffer_size, 1, MPI_INT, rr, 1, MPI_COMM_WORLD, MPI_STATUS_IGNORE); + // while (chunk_buffer_size != 0) { + // char *chunk_buffer = new char[chunk_buffer_size]; + // MPI_Recv(chunk_buffer, chunk_buffer_size, MPI_CHAR, rr, 0, MPI_COMM_WORLD, MPI_STATUS_IGNORE); + // // fputs(chunk_buffer, output); + // int index = 0, last_char = 0; + // char c; + // while (last_char < chunk_buffer_size) { + // c = chunk_buffer[last_char++]; + // virtual_line[index++] = c; + // if (c == '\n') { + // virtual_line[index] = '\0'; + // index = 0; + // p_output->append_line(virtual_line); + // } + // } + // delete[] chunk_buffer; + // MPI_Recv(&chunk_buffer_size, 1, MPI_INT, rr, 1, MPI_COMM_WORLD, MPI_STATUS_IGNORE); + // } + // // if (ri<remotethreads-1) fprintf(output, "\n"); - // now get the binary local file - long buffer_size = 0; - // get the size of the buffer - MPI_Recv(&buffer_size, 1, MPI_LONG, rr, 1, MPI_COMM_WORLD, MPI_STATUS_IGNORE); - // allocate the bufer - char *binary_buffer = new char[buffer_size]; - // actually receive the buffer - MPI_Recv(binary_buffer, buffer_size, MPI_CHAR, rr, 0, MPI_COMM_WORLD, MPI_STATUS_IGNORE); - // we can write it to disk - tppoan.write(binary_buffer, buffer_size); - delete[] binary_buffer; - } + // // now get the binary local file + // long buffer_size = 0; + // // get the size of the buffer + // MPI_Recv(&buffer_size, 1, MPI_LONG, rr, 1, MPI_COMM_WORLD, MPI_STATUS_IGNORE); + // // allocate the bufer + // char *binary_buffer = new char[buffer_size]; + // // actually receive the buffer + // MPI_Recv(binary_buffer, buffer_size, MPI_CHAR, rr, 0, MPI_COMM_WORLD, MPI_STATUS_IGNORE); + // // we can write it to disk + // // tppoan.write(binary_buffer, buffer_size); + // delete[] binary_buffer; + // } } } + #endif #ifdef USE_NVTX nvtxRangePop(); #endif - tppoanp->close(); - delete tppoanp; - } else { // In case TPPOAN could not be opened. Should never happen. - logger->err("\nERROR: failed to open TPPOAN file.\n"); - } - fclose(output); + // tppoanp->close(); + // delete tppoanp; + // } else { // In case TPPOAN could not be opened. Should never happen. + // logger->err("\nERROR: failed to open TPPOAN file.\n"); + // } + // fclose(output); + // p_output->write_to_disk(output_path + "/c_OCLU"); + // delete p_output; // Clean memory delete cid; delete p_scattering_angles; @@ -448,6 +557,8 @@ void cluster(const string& config_file, const string& data_file, const string& o // tppoan.open(tppoan_name.c_str(), ios::out | ios::binary); // Create this variable and initialise it with a default here, so that it is defined anyway, with or without OpenMP support enabled int ompnumthreads = 1; + VirtualAsciiFile **p_outarray = NULL; + VirtualBinaryFile **vtppoanarray = NULL; #pragma omp parallel { @@ -458,10 +569,16 @@ void cluster(const string& config_file, const string& data_file, const string& o myompthread = omp_get_thread_num(); if (myompthread == 0) ompnumthreads = omp_get_num_threads(); #endif + if (myompthread == 0) { + p_outarray = new VirtualAsciiFile*[ompnumthreads]; + vtppoanarray = new VirtualBinaryFile*[ompnumthreads]; + } // To test parallelism, I will now start feeding this function with "clean" copies of the parameters, so that they will not be changed by previous iterations, and each one will behave as the first one. Define all (empty) variables here, so they have the correct scope, then they get different definitions depending on thread number ClusterIterationData *cid_2 = NULL; - FILE *output_2 = NULL; - fstream *tppoanp_2 = NULL; + //FILE *output_2 = NULL; + VirtualAsciiFile *p_output_2 = NULL; + VirtualBinaryFile *vtppoanp_2 = NULL; + // fstream *tppoanp_2 = NULL; // for threads other than the 0, create distinct copies of all relevant data, while for thread 0 just define new references / pointers to the original ones if (myompthread == 0) { cid_2 = cid; @@ -471,17 +588,21 @@ void cluster(const string& config_file, const string& data_file, const string& o // this is not thread 0, so do create fresh copies of all local variables cid_2 = new ClusterIterationData(*cid); } - output_2 = fopen((output_path + "/c_OCLU_" + to_string(mpidata->rank) + "_" + to_string(myompthread)).c_str(), "w"); - tppoanp_2 = new fstream; - tppoanp_2->open((output_path + "/c_TPPOAN_" + to_string(mpidata->rank) + "_" + to_string(myompthread)).c_str(), ios::out | ios::binary); - fstream &tppoan_2 = *tppoanp_2; + // output_2 = fopen((output_path + "/c_OCLU_" + to_string(mpidata->rank) + "_" + to_string(myompthread)).c_str(), "w"); + p_output_2 = new VirtualAsciiFile(); + p_outarray[myompthread] = p_output_2; + vtppoanp_2 = new VirtualBinaryFile(); + vtppoanarray[myompthread] = vtppoanp_2; + // tppoanp_2 = new fstream; + // tppoanp_2->open((output_path + "/c_TPPOAN_" + to_string(mpidata->rank) + "_" + to_string(myompthread)).c_str(), ios::out | ios::binary); + // fstream &tppoan_2 = *tppoanp_2; // make sure all threads align here: I don't want the following loop to accidentally start for thread 0, possibly modifying some variables before they are copied by all other threads #pragma omp barrier if (myompthread==0) logger->log("Syncing OpenMP threads and starting the loop on wavelengths\n"); // ok, now I can actually start the parallel calculations #pragma omp for for (int jxi488 = cid_2->firstxi; jxi488 <= cid_2->lastxi; jxi488++) { - int jer = cluster_jxi488_cycle(jxi488, sconf, gconf, p_scattering_angles, cid_2, output_2, output_path, *tppoanp_2); + int jer = cluster_jxi488_cycle(jxi488, sconf, gconf, p_scattering_angles, cid_2, p_output_2, output_path, vtppoanp_2); } #pragma omp barrier @@ -489,9 +610,11 @@ void cluster(const string& config_file, const string& data_file, const string& o if (myompthread != 0) { delete cid_2; } - fclose(output_2); - tppoanp_2->close(); - delete tppoanp_2; + // fclose(output_2); + // p_output_2->write_to_disk(output_path + "/c_OCLU_" + to_string(mpidata->rank) + "_" + to_string(myompthread)); + // delete p_output_2; + // tppoanp_2->close(); + // delete tppoanp_2; #pragma omp barrier { string message = "INFO: Closing thread-local output files of thread " + to_string(myompthread) + " and syncing threads.\n"; @@ -500,75 +623,87 @@ void cluster(const string& config_file, const string& data_file, const string& o } // closes pragma omp parallel #pragma omp barrier { - // tell MPI process 0 how many threads we have on this process (not necessarily the same across all processes) - MPI_Send(&ompnumthreads, 1, MPI_INT, 0, 1, MPI_COMM_WORLD); - // reopen local files, send them all to MPI process 0 - for (int ri = 0; ri < ompnumthreads; ri++) { - string partial_file_name = output_path + "/c_OCLU_" + to_string(mpidata->rank) + "_" + to_string(ri); - string message = "Copying ASCII output in MPI process " + to_string(mpidata->rank) + " of thread " + to_string(ri) + " of " + to_string(ompnumthreads - 1) + "... "; - logger->log(message, LOG_DEBG); - fstream partial_output; - partial_output.open(partial_file_name.c_str(), ios::in); - partial_output.seekg(0, ios::end); - const long partial_output_size = partial_output.tellg(); - partial_output.close(); - partial_output.open(partial_file_name.c_str(), ios::in); - int chunk_buffer_size = 25165824; // Length of char array with 24Mb size - char *chunk_buffer = new char[chunk_buffer_size](); - int full_chunks = (int)(partial_output_size / chunk_buffer_size); - for (int fi = 0; fi < full_chunks; fi++) { - partial_output.read(chunk_buffer, chunk_buffer_size); - // If EOF is reached, do not send EOF character. - long ptr_position = partial_output.tellg(); - if (ptr_position == partial_output_size) { - chunk_buffer[chunk_buffer_size - 1] = '\0'; - } - // Send the size of the buffer that is being transmitted (Node-0 does not know whether it is full or not) - MPI_Send(&chunk_buffer_size, 1, MPI_INT, 0, 1, MPI_COMM_WORLD); - // Actually send the file contents to Node-0 - MPI_Send(chunk_buffer, chunk_buffer_size, MPI_CHAR, 0, 0, MPI_COMM_WORLD); - } - long ptr_position = partial_output.tellg(); - if (ptr_position < partial_output_size) { - // Send the last partial buffer - chunk_buffer_size = partial_output_size - ptr_position; - delete[] chunk_buffer; - chunk_buffer = new char[chunk_buffer_size]; - partial_output.read(chunk_buffer, chunk_buffer_size); - // chunk_buffer[chunk_buffer_size - 1] = '\0'; - // Send the size of the buffer that is being transmitted (Node-0 does not know whether it is full or not) - MPI_Send(&chunk_buffer_size, 1, MPI_INT, 0, 1, MPI_COMM_WORLD); - // Actually send the file contents to Node-0 - MPI_Send(chunk_buffer, chunk_buffer_size, MPI_CHAR, 0, 0, MPI_COMM_WORLD); - } - // Send a size 0 flag to inform Node-0 that the transmission is complete - chunk_buffer_size = 0; - MPI_Send(&chunk_buffer_size, 1, MPI_INT, 0, 1, MPI_COMM_WORLD); - partial_output.close(); - delete[] chunk_buffer; - remove(partial_file_name.c_str()); - logger->log("done.\n", LOG_DEBG); - - partial_file_name = output_path + "/c_TPPOAN_" + to_string(mpidata->rank) + "_" + to_string(ri); - message = "Copying binary output in MPI process " + to_string(mpidata->rank) + " of thread " + to_string(ri) + " of " + to_string(ompnumthreads - 1) + "... "; - logger->log(message, LOG_DEBG); - fstream partial_tppoan; - partial_tppoan.open(partial_file_name.c_str(), ios::in | ios::binary); - partial_tppoan.seekg(0, ios::end); - long buffer_size = partial_tppoan.tellg(); - char *binary_buffer = new char[buffer_size]; - partial_tppoan.seekg(0, ios::beg); - partial_tppoan.read(binary_buffer, buffer_size); - // tell MPI process 0 how large is the buffer - MPI_Send(&buffer_size, 1, MPI_LONG, 0, 1, MPI_COMM_WORLD); - // actually send the buffer - MPI_Send(binary_buffer, buffer_size, MPI_CHAR, 0, 0, MPI_COMM_WORLD); - // tppoan.write(binary_buffer, buffer_size); - partial_tppoan.close(); - delete[] binary_buffer; - remove(partial_file_name.c_str()); - logger->log("done.\n", LOG_DEBG); + for (int ti=1; ti<ompnumthreads; ti++) { + p_outarray[0]->append(*(p_outarray[ti])); + delete p_outarray[ti]; + vtppoanarray[0]->append(*(vtppoanarray[ti])); + delete vtppoanarray[ti]; } + p_outarray[0]->mpisend(mpidata); + delete p_outarray[0]; + delete[] p_outarray; + vtppoanarray[0]->mpisend(mpidata); + delete vtppoanarray[0]; + delete[] vtppoanarray; + // // tell MPI process 0 how many threads we have on this process (not necessarily the same across all processes) + // MPI_Send(&ompnumthreads, 1, MPI_INT, 0, 1, MPI_COMM_WORLD); + // // reopen local files, send them all to MPI process 0 + // for (int ri = 0; ri < ompnumthreads; ri++) { + // string partial_file_name = output_path + "/c_OCLU_" + to_string(mpidata->rank) + "_" + to_string(ri); + // string message = "Copying ASCII output in MPI process " + to_string(mpidata->rank) + " of thread " + to_string(ri) + " of " + to_string(ompnumthreads - 1) + "... "; + // logger->log(message, LOG_DEBG); + // fstream partial_output; + // partial_output.open(partial_file_name.c_str(), ios::in); + // partial_output.seekg(0, ios::end); + // const long partial_output_size = partial_output.tellg(); + // partial_output.close(); + // partial_output.open(partial_file_name.c_str(), ios::in); + // int chunk_buffer_size = 25165824; // Length of char array with 24Mb size + // char *chunk_buffer = new char[chunk_buffer_size](); + // int full_chunks = (int)(partial_output_size / chunk_buffer_size); + // for (int fi = 0; fi < full_chunks; fi++) { + // partial_output.read(chunk_buffer, chunk_buffer_size); + // // If EOF is reached, do not send EOF character. + // long ptr_position = partial_output.tellg(); + // if (ptr_position == partial_output_size) { + // chunk_buffer[chunk_buffer_size - 1] = '\0'; + // } + // // Send the size of the buffer that is being transmitted (Node-0 does not know whether it is full or not) + // MPI_Send(&chunk_buffer_size, 1, MPI_INT, 0, 1, MPI_COMM_WORLD); + // // Actually send the file contents to Node-0 + // MPI_Send(chunk_buffer, chunk_buffer_size, MPI_CHAR, 0, 0, MPI_COMM_WORLD); + // } + // long ptr_position = partial_output.tellg(); + // if (ptr_position < partial_output_size) { + // // Send the last partial buffer + // chunk_buffer_size = partial_output_size - ptr_position; + // delete[] chunk_buffer; + // chunk_buffer = new char[chunk_buffer_size]; + // partial_output.read(chunk_buffer, chunk_buffer_size); + // // chunk_buffer[chunk_buffer_size - 1] = '\0'; + // // Send the size of the buffer that is being transmitted (Node-0 does not know whether it is full or not) + // MPI_Send(&chunk_buffer_size, 1, MPI_INT, 0, 1, MPI_COMM_WORLD); + // // Actually send the file contents to Node-0 + // MPI_Send(chunk_buffer, chunk_buffer_size, MPI_CHAR, 0, 0, MPI_COMM_WORLD); + // } + // // Send a size 0 flag to inform Node-0 that the transmission is complete + // chunk_buffer_size = 0; + // MPI_Send(&chunk_buffer_size, 1, MPI_INT, 0, 1, MPI_COMM_WORLD); + // partial_output.close(); + // delete[] chunk_buffer; + // remove(partial_file_name.c_str()); + // logger->log("done.\n", LOG_DEBG); + + // string partial_file_name = output_path + "/c_TPPOAN_" + to_string(mpidata->rank) + "_" + to_string(ri); + // string message = "Copying binary output in MPI process " + to_string(mpidata->rank) + " of thread " + to_string(ri) + " of " + to_string(ompnumthreads - 1) + "... "; + // logger->log(message, LOG_DEBG); + // fstream partial_tppoan; + // partial_tppoan.open(partial_file_name.c_str(), ios::in | ios::binary); + // partial_tppoan.seekg(0, ios::end); + // long buffer_size = partial_tppoan.tellg(); + // char *binary_buffer = new char[buffer_size]; + // partial_tppoan.seekg(0, ios::beg); + // partial_tppoan.read(binary_buffer, buffer_size); + // // tell MPI process 0 how large is the buffer + // MPI_Send(&buffer_size, 1, MPI_LONG, 0, 1, MPI_COMM_WORLD); + // // actually send the buffer + // MPI_Send(binary_buffer, buffer_size, MPI_CHAR, 0, 0, MPI_COMM_WORLD); + // // // tppoan.write(binary_buffer, buffer_size); + // partial_tppoan.close(); + // delete[] binary_buffer; + // remove(partial_file_name.c_str()); + // logger->log("done.\n", LOG_DEBG); + // } } // Clean memory delete cid; @@ -587,9 +722,10 @@ void cluster(const string& config_file, const string& data_file, const string& o delete logger; } -int cluster_jxi488_cycle(int jxi488, ScattererConfiguration *sconf, GeometryConfiguration *gconf, ScatteringAngles *sa, ClusterIterationData *cid, FILE *output, const string& output_path, fstream& tppoan) +int cluster_jxi488_cycle(int jxi488, ScattererConfiguration *sconf, GeometryConfiguration *gconf, ScatteringAngles *sa, ClusterIterationData *cid, VirtualAsciiFile *output, const string& output_path, VirtualBinaryFile *vtppoanp) { int nxi = sconf->number_of_scales; + char virtual_line[256]; string message = "INFO: running scale iteration " + to_string(jxi488) + " of " + to_string(nxi) + ".\n"; Logger *logger = new Logger(LOG_DEBG); logger->log(message); @@ -617,7 +753,8 @@ int cluster_jxi488_cycle(int jxi488, ScattererConfiguration *sconf, GeometryConf #ifdef USE_NVTX nvtxRangePush("Prepare matrix calculation"); #endif - fprintf(output, "========== JXI =%3d ====================\n", jxi488); + sprintf(virtual_line, "========== JXI =%3d ====================\n\0", jxi488); + output->append_line(virtual_line); double xi = sconf->get_scale(jxi488 - 1); double exdc = sconf->exdc; double exri = sqrt(exdc); @@ -626,15 +763,18 @@ int cluster_jxi488_cycle(int jxi488, ScattererConfiguration *sconf, GeometryConf if (idfc >= 0) { cid->vk = xi * cid->wn; vkarg = cid->vk; - fprintf(output, " VK=%15.7lE, XI=%15.7lE\n", cid->vk, xi); + sprintf(virtual_line, " VK=%15.7lE, XI=%15.7lE\n\0", cid->vk, xi); + output->append_line(virtual_line); } else { vkarg = xi * cid->vk; cid->sqsfi = 1.0 / (xi * xi); - fprintf(output, " XI=%15.7lE\n", xi); + sprintf(virtual_line, " XI=%15.7lE\n\0", xi); + output->append_line(virtual_line); } hjv(exri, vkarg, jer, lcalc, cid->arg, cid->c1, cid->c1ao, cid->c4); if (jer != 0) { - fprintf(output, " STOP IN HJV\n"); + sprintf(virtual_line, " STOP IN HJV\n\0"); + output->append_line(virtual_line); return jer; // break; // rewrite this to go to the end of the function, to free locally allocated variables and return jer } @@ -663,7 +803,8 @@ int cluster_jxi488_cycle(int jxi488, ScattererConfiguration *sconf, GeometryConf cid->c1, cid->c2, jer, lcalc, cid->arg ); if (jer != 0) { - fprintf(output, " STOP IN DME\n"); + sprintf(virtual_line, " STOP IN DME\n\0"); + output->append_line(virtual_line); return jer; //break; } @@ -722,9 +863,11 @@ int cluster_jxi488_cycle(int jxi488, ScattererConfiguration *sconf, GeometryConf } // label 156: continue from here if (inpol == 0) { - fprintf(output, " LIN\n"); + sprintf(virtual_line, " LIN\n\0"); + output->append_line(virtual_line); } else { // label 158 - fprintf(output, " CIRC\n"); + sprintf(virtual_line, " CIRC\n\0"); + output->append_line(virtual_line); } // label 160 double cs0 = 0.25 * cid->vk * cid->vk * cid->vk / acos(0.0); @@ -734,7 +877,10 @@ int cluster_jxi488_cycle(int jxi488, ScattererConfiguration *sconf, GeometryConf double sqk = cid->vk * cid->vk * exdc; aps(cid->zpv, cid->c4->li, nsph, cid->c1, sqk, cid->gaps); rabas(inpol, cid->c4->li, nsph, cid->c1, cid->tqse, cid->tqspe, cid->tqss, cid->tqsps); - if (cid->c4->li != cid->c4->le) fprintf(output, " SPHERES; LMX=LI\n"); + if (cid->c4->li != cid->c4->le) { + sprintf(virtual_line, " SPHERES; LMX=LI\n\0"); + output->append_line(virtual_line); + } for (int i170 = 1; i170 <= nsph; i170++) { if (cid->c1->iog[i170 - 1] >= i170) { int i = i170 - 1; @@ -742,39 +888,57 @@ int cluster_jxi488_cycle(int jxi488, ScattererConfiguration *sconf, GeometryConf cid->c1->sqscs[i] *= cid->sqsfi; cid->c1->sqabs[i] *= cid->sqsfi; cid->c1->sqexs[i] *= cid->sqsfi; - fprintf(output, " SPHERE %2d\n", i170); + sprintf(virtual_line, " SPHERE %2d\n\0", i170); + output->append_line(virtual_line); if (cid->c1->nshl[i] != 1) { - fprintf(output, " SIZE=%15.7lE\n", cid->c2->vsz[i]); + sprintf(virtual_line, " SIZE=%15.7lE\n\0", cid->c2->vsz[i]); + output->append_line(virtual_line); } else { // label 162 - fprintf(output, " SIZE=%15.7lE, REFRACTIVE INDEX=%15.7lE%15.7lE\n", cid->c2->vsz[i], real(cid->c2->vkt[i]), imag(cid->c2->vkt[i])); + sprintf(virtual_line, " SIZE=%15.7lE, REFRACTIVE INDEX=%15.7lE%15.7lE\n\0", cid->c2->vsz[i], real(cid->c2->vkt[i]), imag(cid->c2->vkt[i])); + output->append_line(virtual_line); } // label 164 - fprintf(output, " ----- SCS ----- ABS ----- EXS ----- ALBEDS --\n"); - fprintf(output, " %14.7lE%15.7lE%15.7lE%15.7lE\n", cid->c1->sscs[i], cid->c1->sabs[i], cid->c1->sexs[i], albeds); - fprintf(output, " ---- SCS/GS -- ABS/GS -- EXS/GS ---\n"); - fprintf(output, " %14.7lE%15.7lE%15.7lE\n", cid->c1->sqscs[i], cid->c1->sqabs[i], cid->c1->sqexs[i]); - fprintf(output, " FSAS=%15.7lE%15.7lE\n", real(cid->c1->fsas[i]), imag(cid->c1->fsas[i])); + sprintf(virtual_line, " ----- SCS ----- ABS ----- EXS ----- ALBEDS --\n\0"); + output->append_line(virtual_line); + sprintf(virtual_line, " %14.7lE%15.7lE%15.7lE%15.7lE\n\0", cid->c1->sscs[i], cid->c1->sabs[i], cid->c1->sexs[i], albeds); + output->append_line(virtual_line); + sprintf(virtual_line, " ---- SCS/GS -- ABS/GS -- EXS/GS ---\n\0"); + output->append_line(virtual_line); + sprintf(virtual_line, " %14.7lE%15.7lE%15.7lE\n\0", cid->c1->sqscs[i], cid->c1->sqabs[i], cid->c1->sqexs[i]); + output->append_line(virtual_line); + sprintf(virtual_line, " FSAS=%15.7lE%15.7lE\n\0", real(cid->c1->fsas[i]), imag(cid->c1->fsas[i])); + output->append_line(virtual_line); + double alamb = 2.0 * 3.141592653589793 / cid->vk; + sprintf(virtual_line, "INSERTION: CS_SPHERE %15.7lE%15.7lE%15.7lE%15.7lE\n\0", alamb, cid->c1->sscs[i], cid->c1->sabs[i], cid->c1->sexs[i]); + output->append_line(virtual_line); csch = 2.0 * cid->vk * cid->sqsfi / cid->c1->gcsv[i]; s0 = cid->c1->fsas[i] * exri; qschu = imag(s0) * csch; pschu = real(s0) * csch; s0mag = cabs(s0) * cs0; - fprintf(output, " QSCHU=%15.7lE, PSCHU=%15.7lE, S0MAG=%15.7lE\n", qschu, pschu, s0mag); + sprintf(virtual_line, " QSCHU=%15.7lE, PSCHU=%15.7lE, S0MAG=%15.7lE\n\0", qschu, pschu, s0mag); + output->append_line(virtual_line); double rapr = cid->c1->sexs[i] - cid->gaps[i]; double cosav = cid->gaps[i] / cid->c1->sscs[i]; - fprintf(output, " COSAV=%15.7lE, RAPRS=%15.7lE\n", cosav, rapr); - fprintf(output, " IPO= 1, TQEk=%15.7lE, TQSk=%15.7lE\n", cid->tqse[0][i], cid->tqss[0][i]); - fprintf(output, " IPO= 2, TQEk=%15.7lE, TQSk=%15.7lE\n", cid->tqse[1][i], cid->tqss[1][i]); + sprintf(virtual_line, " COSAV=%15.7lE, RAPRS=%15.7lE\n\0", cosav, rapr); + output->append_line(virtual_line); + sprintf(virtual_line, " IPO= 1, TQEk=%15.7lE, TQSk=%15.7lE\n\0", cid->tqse[0][i], cid->tqss[0][i]); + output->append_line(virtual_line); + sprintf(virtual_line, " IPO= 2, TQEk=%15.7lE, TQSk=%15.7lE\n\0", cid->tqse[1][i], cid->tqss[1][i]); + output->append_line(virtual_line); } } // i170 loop - fprintf(output, " FSAT=%15.7lE%15.7lE\n", real(cid->c3->tfsas), imag(cid->c3->tfsas)); + sprintf(virtual_line, " FSAT=%15.7lE%15.7lE\n\0", real(cid->c3->tfsas), imag(cid->c3->tfsas)); + output->append_line(virtual_line); csch = 2.0 * cid->vk * cid->sqsfi / cid->c3->gcs; s0 = cid->c3->tfsas * exri; qschu = imag(s0) * csch; pschu = real(s0) * csch; s0mag = cabs(s0) * cs0; - fprintf(output, " QSCHU=%15.7lE, PSCHU=%15.7lE, S0MAG=%15.7lE\n", qschu, pschu, s0mag); - tppoan.write(reinterpret_cast<char *>(&(cid->vk)), sizeof(double)); + sprintf(virtual_line, " QSCHU=%15.7lE, PSCHU=%15.7lE, S0MAG=%15.7lE\n\0", qschu, pschu, s0mag); + output->append_line(virtual_line); + // tppoan.write(reinterpret_cast<char *>(&(cid->vk)), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(cid->vk)); pcrsm0(cid->vk, exri, inpol, cid->c1, cid->c1ao, cid->c4); apcra(cid->zpv, cid->c4->le, cid->c1ao->am0m, inpol, sqk, cid->gapm, cid->gappm); #ifdef USE_NVTX @@ -782,7 +946,7 @@ int cluster_jxi488_cycle(int jxi488, ScattererConfiguration *sconf, GeometryConf #endif interval_end = chrono::high_resolution_clock::now(); elapsed = interval_end - interval_start; - message = "INFO: average calculation for scale " + to_string(jxi488) + " took " + to_string(elapsed.count()) + "s.\n"; + message = "INFO: average calculation for scale " + to_string(jxi488) + " took " + to_string(elapsed.count()) + "s.\n\0"; logger->log(message); interval_start = chrono::high_resolution_clock::now(); #ifdef USE_NVTX @@ -875,11 +1039,16 @@ int cluster_jxi488_cycle(int jxi488, ScattererConfiguration *sconf, GeometryConf jw = 1; } // label 196 - tppoan.write(reinterpret_cast<char *>(&th), sizeof(double)); - tppoan.write(reinterpret_cast<char *>(&ph), sizeof(double)); - tppoan.write(reinterpret_cast<char *>(&ths), sizeof(double)); - tppoan.write(reinterpret_cast<char *>(&phs), sizeof(double)); - tppoan.write(reinterpret_cast<char *>(&(cid->scan)), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&th), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(th)); + // tppoan.write(reinterpret_cast<char *>(&ph), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(ph)); + // tppoan.write(reinterpret_cast<char *>(&ths), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(ths)); + // tppoan.write(reinterpret_cast<char *>(&phs), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(phs)); + // tppoan.write(reinterpret_cast<char *>(&(cid->scan)), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(cid->scan)); if (jaw != 0) { jaw = 0; mextc(cid->vk, exri, cid->c1ao->fsacm, cid->cextlr, cid->cext); @@ -887,34 +1056,45 @@ int cluster_jxi488_cycle(int jxi488, ScattererConfiguration *sconf, GeometryConf for (int i = 0; i < 4; i++) { for (int j = 0; j < 4; j++) { double value = cid->cext[i][j]; - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); } } for (int i = 0; i < 2; i++) { double value = cid->c1ao->scscm[i]; - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); value = real(cid->c1ao->scscpm[i]); - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); value = imag(cid->c1ao->scscpm[i]); - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); value = cid->c1ao->ecscm[i]; - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); value = real(cid->c1ao->ecscpm[i]); - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); value = imag(cid->c1ao->ecscpm[i]); - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); } for (int i = 0; i < 3; i++) { for (int j = 0; j < 2; j++) { double value = cid->gapm[i][j]; - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); value = real(cid->gappm[i][j]); - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); value = imag(cid->gappm[i][j]); - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); } } - fprintf(output, " CLUSTER (ENSEMBLE AVERAGE, MODE%2d)\n", iavm); + sprintf(virtual_line, " CLUSTER (ENSEMBLE AVERAGE, MODE%2d)\n\0", iavm); + output->append_line(virtual_line); int jlr = 2; for (int ilr210 = 1; ilr210 <= 2; ilr210++) { int ipol = (ilr210 % 2 == 0) ? 1 : -1; @@ -938,104 +1118,145 @@ int cluster_jxi488_cycle(int jxi488, ScattererConfiguration *sconf, GeometryConf double rfinrm = real(cid->c1ao->fsacm[ilr210 - 1][ilr210 - 1]) / real(cid->c3->tfsas); double extcrm = imag(cid->c1ao->fsacm[ilr210 - 1][ilr210 - 1]) / imag(cid->c3->tfsas); if (inpol == 0) { - fprintf(output, " LIN %2d\n", ipol); + sprintf(virtual_line, " LIN %2d\n\0", ipol); + output->append_line(virtual_line); } else { // label 206 - fprintf(output, " CIRC %2d\n", ipol); + sprintf(virtual_line, " CIRC %2d\n\0", ipol); + output->append_line(virtual_line); } // label 208 - fprintf(output, " ----- SCC ----- ABC ----- EXC ----- ALBEDC --\n"); - fprintf(output, " %14.7lE%15.7lE%15.7lE%15.7lE\n", scasm, abssm, extsm, albdm); - fprintf(output, " --- SCC/TGS - ABC/TGS - EXC/TGS ---\n"); - fprintf(output, " %14.7lE%15.7lE%15.7lE\n", qscam, qabsm, qextm); - fprintf(output, " ---- SCCRT --- ABCRT --- EXCRT ----\n"); - fprintf(output, " %14.7lE%15.7lE%15.7lE\n", scarm, absrm, extrm); - fprintf( - output, " FSAC(%1d,%1d)=%15.7lE%15.7lE FSAC(%1d,%1d)=%15.7lE%15.7lE\n", + sprintf(virtual_line, " ----- SCC ----- ABC ----- EXC ----- ALBEDC --\n\0"); + output->append_line(virtual_line); + sprintf(virtual_line, " %14.7lE%15.7lE%15.7lE%15.7lE\n\0", scasm, abssm, extsm, albdm); + output->append_line(virtual_line); + sprintf(virtual_line, " --- SCC/TGS - ABC/TGS - EXC/TGS ---\n\0"); + output->append_line(virtual_line); + sprintf(virtual_line, " %14.7lE%15.7lE%15.7lE\n\0", qscam, qabsm, qextm); + output->append_line(virtual_line); + sprintf(virtual_line, " ---- SCCRT --- ABCRT --- EXCRT ----\n\0"); + output->append_line(virtual_line); + sprintf(virtual_line, " %14.7lE%15.7lE%15.7lE\n\0", scarm, absrm, extrm); + output->append_line(virtual_line); + sprintf( + virtual_line, " FSAC(%1d,%1d)=%15.7lE%15.7lE FSAC(%1d,%1d)=%15.7lE%15.7lE\n\0", ilr210, ilr210, real(cid->c1ao->fsacm[ilr210 - 1][ilr210 - 1]), imag(cid->c1ao->fsacm[ilr210 - 1][ilr210 - 1]), jlr, ilr210, real(cid->c1ao->fsacm[jlr - 1][ilr210 - 1]), imag(cid->c1ao->fsacm[jlr - 1][ilr210 - 1]) ); - fprintf( - output, " RE(FSAC(%1d,%1d))/RE(TFSAS)=%15.7lE, IM(FSAC(%1d,%1d))/IM(TFSAS)=%15.7lE\n", + output->append_line(virtual_line); + sprintf( + virtual_line, " RE(FSAC(%1d,%1d))/RE(TFSAS)=%15.7lE, IM(FSAC(%1d,%1d))/IM(TFSAS)=%15.7lE\n\0", ilr210, ilr210, rfinrm, ilr210, ilr210, extcrm ); - fprintf(output, " QSCHU=%15.7lE, PSCHU=%15.7lE, S0MAG=%15.7lE\n", qschum, pschum, s0magm); + output->append_line(virtual_line); + sprintf(virtual_line, " QSCHU=%15.7lE, PSCHU=%15.7lE, S0MAG=%15.7lE\n\0", qschum, pschum, s0magm); + output->append_line(virtual_line); double rapr = cid->c1ao->ecscm[ilr210 - 1] - cid->gapm[2][ilr210 - 1]; double cosav = cid->gapm[2][ilr210 - 1] / cid->c1ao->scscm[ilr210 - 1]; double fz = rapr; - fprintf(output, " COSAV=%15.7lE, RAPRS=%15.7lE\n", cosav, rapr); - fprintf(output, " Fk=%15.7lE\n", fz); + sprintf(virtual_line, " COSAV=%15.7lE, RAPRS=%15.7lE\n\0", cosav, rapr); + output->append_line(virtual_line); + sprintf(virtual_line, " Fk=%15.7lE\n\0", fz); + output->append_line(virtual_line); + double alamb = 2.0 * 3.141592653589793 / cid->vk; + sprintf(virtual_line, "INSERTION: CSM_CLUSTER %15.7lE%15.7lE%15.7lE%15.7lE\n\0", alamb, scasm, abssm, extsm); + output->append_line(virtual_line); } // ilr210 loop double rmbrif = (real(cid->c1ao->fsacm[0][0]) - real(cid->c1ao->fsacm[1][1])) / real(cid->c1ao->fsacm[0][0]); double rmdchr = (imag(cid->c1ao->fsacm[0][0]) - imag(cid->c1ao->fsacm[1][1])) / imag(cid->c1ao->fsacm[0][0]); - fprintf(output, " (RE(FSAC(1,1))-RE(FSAC(2,2)))/RE(FSAC(1,1))=%15.7lE\n", rmbrif); - fprintf(output, " (IM(FSAC(1,1))-IM(FSAC(2,2)))/IM(FSAC(1,1))=%15.7lE\n", rmdchr); + sprintf(virtual_line, " (RE(FSAC(1,1))-RE(FSAC(2,2)))/RE(FSAC(1,1))=%15.7lE\n\0", rmbrif); + output->append_line(virtual_line); + sprintf(virtual_line, " (IM(FSAC(1,1))-IM(FSAC(2,2)))/IM(FSAC(1,1))=%15.7lE\n\0", rmdchr); + output->append_line(virtual_line); } // label 212 - fprintf(output, "********** JTH =%3d, JPH =%3d, JTHS =%3d, JPHS =%3d ********************\n", jth486, jph484, jths, jphs); - fprintf(output, " TIDG=%10.3lE, PIDG=%10.3lE, TSDG=%10.3lE, PSDG=%10.3lE\n", th, ph, ths, phs); - fprintf(output, " SCAND=%10.3lE\n", cid->scan); - fprintf(output, " CFMP=%15.7lE, SFMP=%15.7lE\n", cid->cfmp, cid->sfmp); - fprintf(output, " CFSP=%15.7lE, SFSP=%15.7lE\n", cid->cfsp, cid->sfsp); + sprintf(virtual_line, "********** JTH =%3d, JPH =%3d, JTHS =%3d, JPHS =%3d ********************\n\0", jth486, jph484, jths, jphs); + output->append_line(virtual_line); + sprintf(virtual_line, " TIDG=%10.3lE, PIDG=%10.3lE, TSDG=%10.3lE, PSDG=%10.3lE\n\0", th, ph, ths, phs); + output->append_line(virtual_line); + sprintf(virtual_line, " SCAND=%10.3lE\n\0", cid->scan); + output->append_line(virtual_line); + sprintf(virtual_line, " CFMP=%15.7lE, SFMP=%15.7lE\n\0", cid->cfmp, cid->sfmp); + output->append_line(virtual_line); + sprintf(virtual_line, " CFSP=%15.7lE, SFSP=%15.7lE\n\0", cid->cfsp, cid->sfsp); + output->append_line(virtual_line); if (isam >= 0) { - fprintf(output, " UNI=(%12.5lE,%12.5lE,%12.5lE)\n", cid->un[0], cid->un[1], cid->un[2]); - fprintf(output, " UNS=(%12.5lE,%12.5lE,%12.5lE)\n", cid->uns[0], cid->uns[1], cid->uns[2]); + sprintf(virtual_line, " UNI=(%12.5lE,%12.5lE,%12.5lE)\n\0", cid->un[0], cid->un[1], cid->un[2]); + output->append_line(virtual_line); + sprintf(virtual_line, " UNS=(%12.5lE,%12.5lE,%12.5lE)\n\0", cid->uns[0], cid->uns[1], cid->uns[2]); + output->append_line(virtual_line); } else { // label 214 - fprintf(output, " UN=(%12.5lE,%12.5lE,%12.5lE)\n\n", cid->un[0], cid->un[1], cid->un[2]); + sprintf(virtual_line, " UN=(%12.5lE,%12.5lE,%12.5lE)\n\n\0", cid->un[0], cid->un[1], cid->un[2]); + output->append_line(virtual_line); } // label 220 if (inpol == 0) { - fprintf(output, " LIN\n"); + sprintf(virtual_line, " LIN\n\0"); + output->append_line(virtual_line); } else { // label 222 - fprintf(output, " CIRC\n"); + sprintf(virtual_line, " CIRC\n\0"); + output->append_line(virtual_line); } // label 224 scr2(cid->vk, vkarg, exri, cid->duk, cid->c1, cid->c1ao, cid->c3, cid->c4); - if (cid->c4->li != cid->c4->le) fprintf(output, " SPHERES; LMX=MIN0(LI,LE)\n"); + if (cid->c4->li != cid->c4->le) { + sprintf(virtual_line, " SPHERES; LMX=MIN0(LI,LE)\n\0"); + output->append_line(virtual_line); + } for (int i226 = 1; i226 <= nsph; i226++) { if (cid->c1->iog[i226 - 1] >= i226) { - fprintf(output, " SPHERE %2d\n", i226); - fprintf( - output, " SAS(1,1)=%15.7lE%15.7lE, SAS(2,1)=%15.7lE%15.7lE\n", + sprintf(virtual_line, " SPHERE %2d\n\0", i226); + output->append_line(virtual_line); + sprintf( + virtual_line, " SAS(1,1)=%15.7lE%15.7lE, SAS(2,1)=%15.7lE%15.7lE\n\0", real(cid->c1->sas[i226 - 1][0][0]), imag(cid->c1->sas[i226 - 1][0][0]), real(cid->c1->sas[i226 - 1][1][0]), imag(cid->c1->sas[i226 - 1][1][0]) ); - fprintf( - output, " SAS(1,2)=%15.7lE%15.7lE, SAS(2,2)=%15.7lE%15.7lE\n", + output->append_line(virtual_line); + sprintf( + virtual_line, " SAS(1,2)=%15.7lE%15.7lE, SAS(2,2)=%15.7lE%15.7lE\n\0", real(cid->c1->sas[i226 - 1][0][1]), imag(cid->c1->sas[i226 - 1][0][1]), real(cid->c1->sas[i226 - 1][1][1]), imag(cid->c1->sas[i226 - 1][1][1]) ); + output->append_line(virtual_line); for (int j225 = 0; j225 < 16; j225++) { cid->c1->vint[j225] = cid->c1->vints[i226 - 1][j225]; } // j225 loop mmulc(cid->c1->vint, cid->cmullr, cid->cmul); - fprintf(output, " MULS\n"); + sprintf(virtual_line, " MULS\n\0"); + output->append_line(virtual_line); for (int i1 = 0; i1 < 4; i1++) { - fprintf( - output, " %15.7lE%15.7lE%15.7lE%15.7lE\n", + sprintf( + virtual_line, " %15.7lE%15.7lE%15.7lE%15.7lE\n\0", cid->cmul[i1][0], cid->cmul[i1][1], cid->cmul[i1][2], cid->cmul[i1][3] ); + output->append_line(virtual_line); } // i1 loop - fprintf(output, " MULSLR\n"); + sprintf(virtual_line, " MULSLR\n\0"); + output->append_line(virtual_line); for (int i1 = 0; i1 < 4; i1++) { - fprintf( - output, " %15.7lE%15.7lE%15.7lE%15.7lE\n", + sprintf( + virtual_line, " %15.7lE%15.7lE%15.7lE%15.7lE\n\0", cid->cmullr[i1][0], cid->cmullr[i1][1], cid->cmullr[i1][2], cid->cmullr[i1][3] ); + output->append_line(virtual_line); } // i1 loop } } // i226 loop - fprintf( - output, " SAT(1,1)=%15.7lE%15.7lE, SAT(2,1)=%15.7lE%15.7lE\n", + sprintf( + virtual_line, " SAT(1,1)=%15.7lE%15.7lE, SAT(2,1)=%15.7lE%15.7lE\n\0", real(cid->c3->tsas[0][0]), imag(cid->c3->tsas[0][0]), real(cid->c3->tsas[1][0]), imag(cid->c3->tsas[1][0]) ); - fprintf( - output, " SAT(1,2)=%15.7lE%15.7lE, SAT(2,2)=%15.7lE%15.7lE\n", + output->append_line(virtual_line); + sprintf( + virtual_line, " SAT(1,2)=%15.7lE%15.7lE, SAT(2,2)=%15.7lE%15.7lE\n\0", real(cid->c3->tsas[0][1]), imag(cid->c3->tsas[0][1]), real(cid->c3->tsas[1][1]), imag(cid->c3->tsas[1][1]) ); - fprintf(output, " CLUSTER\n"); + output->append_line(virtual_line); + sprintf(virtual_line, " CLUSTER\n\0"); + output->append_line(virtual_line); pcros(cid->vk, exri, cid->c1, cid->c1ao, cid->c4); mextc(cid->vk, exri, cid->c1ao->fsac, cid->cextlr, cid->cext); mmulc(cid->c1->vint, cid->cmullr, cid->cmul); @@ -1045,73 +1266,95 @@ int cluster_jxi488_cycle(int jxi488, ScattererConfiguration *sconf, GeometryConf for (int i = 0; i < 4; i++) { for (int j = 0; j < 4; j++) { double value = cid->cext[i][j]; - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); } } for (int i = 0; i < 2; i++) { double value = cid->c1ao->scsc[i]; - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); value = real(cid->c1ao->scscp[i]); - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); value = imag(cid->c1ao->scscp[i]); - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); value = cid->c1ao->ecsc[i]; - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); value = real(cid->c1ao->ecscp[i]); - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); value = imag(cid->c1ao->ecscp[i]); - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); } for (int i = 0; i < 3; i++) { for (int j = 0; j < 2; j++) { double value = cid->gap[i][j]; - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); value = real(cid->gapp[i][j]); - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); value = imag(cid->gapp[i][j]); - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); } } for (int i = 0; i < 2; i++) { for (int j = 0; j < 3; j++) { double value = cid->tqce[i][j]; - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); value = real(cid->tqcpe[i][j]); - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); value = imag(cid->tqcpe[i][j]); - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); } } for (int i = 0; i < 2; i++) { for (int j = 0; j < 3; j++) { double value = cid->tqcs[i][j]; - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); value = real(cid->tqcps[i][j]); - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); value = imag(cid->tqcps[i][j]); - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); } } for (int i = 0; i < 3; i++) { double value = cid->u[i]; - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); value = cid->up[i]; - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); value = cid->un[i]; - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); } } // label 254 for (int i = 0; i < 16; i++) { double value = real(cid->c1->vint[i]); - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); value = imag(cid->c1->vint[i]); - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); } for (int i = 0; i < 4; i++) { for (int j = 0; j < 4; j++) { double value = cid->cmul[i][j]; - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); } } int jlr = 2; @@ -1137,44 +1380,63 @@ int cluster_jxi488_cycle(int jxi488, ScattererConfiguration *sconf, GeometryConf double refinr = real(cid->c1ao->fsac[ilr290 - 1][ilr290 - 1]) / real(cid->c3->tfsas); double extcor = imag(cid->c1ao->fsac[ilr290 - 1][ilr290 - 1]) / imag(cid->c3->tfsas); if (inpol == 0) { - fprintf(output, " LIN %2d\n", ipol); + sprintf(virtual_line, " LIN %2d\n\0", ipol); + output->append_line(virtual_line); } else { // label 273 - fprintf(output, " CIRC %2d\n", ipol); + sprintf(virtual_line, " CIRC %2d\n\0", ipol); + output->append_line(virtual_line); } // label 275 - fprintf(output, " ----- SCC ----- ABC ----- EXC ----- ALBEDC --\n"); - fprintf( - output, " %14.7lE%15.7lE%15.7lE%15.7lE\n", + sprintf(virtual_line, " ----- SCC ----- ABC ----- EXC ----- ALBEDC --\n\0"); + output->append_line(virtual_line); + sprintf( + virtual_line, " %14.7lE%15.7lE%15.7lE%15.7lE\n\0", scasec, abssec, extsec, albedc ); - fprintf(output, " --- SCC/TGS - ABC/TGS - EXC/TGS ---\n"); - fprintf( - output, " %14.7lE%15.7lE%15.7lE\n", + output->append_line(virtual_line); + sprintf(virtual_line, " --- SCC/TGS - ABC/TGS - EXC/TGS ---\n\0"); + output->append_line(virtual_line); + sprintf( + virtual_line, " %14.7lE%15.7lE%15.7lE\n\0", qsca, qabs, qext ); - fprintf(output, " ---- SCCRT --- ABCRT --- EXCRT ----\n"); - fprintf( - output, " %14.7lE%15.7lE%15.7lE\n", + output->append_line(virtual_line); + sprintf(virtual_line, " ---- SCCRT --- ABCRT --- EXCRT ----\n\0"); + output->append_line(virtual_line); + sprintf( + virtual_line, " %14.7lE%15.7lE%15.7lE\n\0", scarat, absrat, extrat ); - fprintf( - output, " FSAC(%1d,%1d)=%15.7lE%15.7lE FSAC(%1d,%1d)=%15.7lE%15.7lE\n", + output->append_line(virtual_line); + sprintf( + virtual_line, " FSAC(%1d,%1d)=%15.7lE%15.7lE FSAC(%1d,%1d)=%15.7lE%15.7lE\n\0", ilr290, ilr290, real(cid->c1ao->fsac[ilr290 - 1][ilr290 - 1]), imag(cid->c1ao->fsac[ilr290 - 1][ilr290 - 1]), jlr, ilr290, real(cid->c1ao->fsac[jlr - 1][ilr290 - 1]), imag(cid->c1ao->fsac[jlr - 1][ilr290 - 1]) ); - fprintf( - output, " SAC(%1d,%1d)=%15.7lE%15.7lE SAC(%1d,%1d)=%15.7lE%15.7lE\n", + output->append_line(virtual_line); + sprintf( + virtual_line, " SAC(%1d,%1d)=%15.7lE%15.7lE SAC(%1d,%1d)=%15.7lE%15.7lE\n\0", ilr290, ilr290, real(cid->c1ao->sac[ilr290 - 1][ilr290 - 1]), imag(cid->c1ao->sac[ilr290 - 1][ilr290 - 1]), jlr, ilr290, real(cid->c1ao->sac[jlr - 1][ilr290 - 1]), imag(cid->c1ao->sac[jlr - 1][ilr290 - 1]) ); - fprintf( - output, " RE(FSAC(%1d,%1d))/RE(TFSAS)=%15.7lE, IM(FSAC(%1d,%1d))/IM(TFSAS)=%15.7lE\n", + output->append_line(virtual_line); + sprintf( + virtual_line, " RE(FSAC(%1d,%1d))/RE(TFSAS)=%15.7lE, IM(FSAC(%1d,%1d))/IM(TFSAS)=%15.7lE\n\0", ilr290, ilr290, refinr, ilr290, ilr290, extcor ); - fprintf( - output, " QSCHU=%15.7lE, PSCHU=%15.7lE, S0MAG=%15.7lE\n", + output->append_line(virtual_line); + sprintf( + virtual_line, " QSCHU=%15.7lE, PSCHU=%15.7lE, S0MAG=%15.7lE\n\0", qschu, pschu, s0mag ); + output->append_line(virtual_line); + double alamb = 2.0 * 3.141592653589793 / cid->vk; + if (ilr290 == 1) { + sprintf(virtual_line, "INSERTION: CS1_CLUSTER %9.2lf%9.2lf%9.2lf%9.2lf%9.2lf%13.5lf%13.5lf%13.5lf\n\0", alamb, th, ph, ths, phs, scasec, abssec, extsec); + } else if (ilr290 == 2) { + sprintf(virtual_line, "INSERTION: CS2_CLUSTER %9.2lf%9.2lf%9.2lf%9.2lf%9.2lf%13.5lf%13.5lf%13.5lf\n\0", alamb, th, ph, ths, phs, scasec, abssec, extsec); + } + output->append_line(virtual_line); bool goto190 = isam >= 0 && (jths > 1 || jphs > 1); if (!goto190) { cid->gapv[0] = cid->gap[0][ilr290 - 1]; @@ -1184,9 +1446,12 @@ int cluster_jxi488_cycle(int jxi488, ScattererConfiguration *sconf, GeometryConf double scatts = cid->c1ao->scsc[ilr290 - 1]; double rapr, cosav, fp, fn, fk, fx, fy, fz; rftr(cid->u, cid->up, cid->un, cid->gapv, extins, scatts, rapr, cosav, fp, fn, fk, fx, fy, fz); - fprintf(output, " COSAV=%15.7lE, RAPRS=%15.7lE\n", cosav, rapr); - fprintf(output, " Fl=%15.7lE, Fr=%15.7lE, Fk=%15.7lE\n", fp, fn, fk); - fprintf(output, " Fx=%15.7lE, Fy=%15.7lE, Fz=%15.7lE\n", fx, fy, fz); + sprintf(virtual_line, " COSAV=%15.7lE, RAPRS=%15.7lE\n\0", cosav, rapr); + output->append_line(virtual_line); + sprintf(virtual_line, " Fl=%15.7lE, Fr=%15.7lE, Fk=%15.7lE\n\0", fp, fn, fk); + output->append_line(virtual_line); + sprintf(virtual_line, " Fx=%15.7lE, Fy=%15.7lE, Fz=%15.7lE\n\0", fx, fy, fz); + output->append_line(virtual_line); cid->tqev[0] = cid->tqce[ilr290 - 1][0]; cid->tqev[1] = cid->tqce[ilr290 - 1][1]; cid->tqev[2] = cid->tqce[ilr290 - 1][2]; @@ -1195,35 +1460,45 @@ int cluster_jxi488_cycle(int jxi488, ScattererConfiguration *sconf, GeometryConf cid->tqsv[2] = cid->tqcs[ilr290 - 1][2]; double tep, ten, tek, tsp, tsn, tsk; tqr(cid->u, cid->up, cid->un, cid->tqev, cid->tqsv, tep, ten, tek, tsp, tsn, tsk); - fprintf(output, " TQEl=%15.7lE, TQEr=%15.7lE, TQEk=%15.7lE\n", tep, ten, tek); - fprintf(output, " TQSl=%15.7lE, TQSr=%15.7lE, TQSk=%15.7lE\n", tsp, tsn, tsk); - fprintf( - output, " TQEx=%15.7lE, TQEy=%15.7lE, TQEz=%15.7lE\n", + sprintf(virtual_line, " TQEl=%15.7lE, TQEr=%15.7lE, TQEk=%15.7lE\n\0", tep, ten, tek); + output->append_line(virtual_line); + sprintf(virtual_line, " TQSl=%15.7lE, TQSr=%15.7lE, TQSk=%15.7lE\n\0", tsp, tsn, tsk); + output->append_line(virtual_line); + sprintf( + virtual_line, " TQEx=%15.7lE, TQEy=%15.7lE, TQEz=%15.7lE\n\0", cid->tqce[ilr290 - 1][0], cid->tqce[ilr290 - 1][1], cid->tqce[ilr290 - 1][2] ); - fprintf( - output, " TQSx=%15.7lE, TQSy=%15.7lE, TQSz=%15.7lE\n", + output->append_line(virtual_line); + sprintf( + virtual_line, " TQSx=%15.7lE, TQSy=%15.7lE, TQSz=%15.7lE\n\0", cid->tqcs[ilr290 - 1][0], cid->tqcs[ilr290 - 1][1], cid->tqcs[ilr290 - 1][2] ); + output->append_line(virtual_line); } } //ilr290 loop double rbirif = (real(cid->c1ao->fsac[0][0]) - real(cid->c1ao->fsac[1][1])) / real(cid->c1ao->fsac[0][0]); double rdichr = (imag(cid->c1ao->fsac[0][0]) - imag(cid->c1ao->fsac[1][1])) / imag(cid->c1ao->fsac[0][0]); - fprintf(output, " (RE(FSAC(1,1))-RE(FSAC(2,2)))/RE(FSAC(1,1))=%15.7lE\n", rbirif); - fprintf(output, " (IM(FSAC(1,1))-IM(FSAC(2,2)))/IM(FSAC(1,1))=%15.7lE\n", rdichr); - fprintf(output, " MULC\n"); + sprintf(virtual_line, " (RE(FSAC(1,1))-RE(FSAC(2,2)))/RE(FSAC(1,1))=%15.7lE\n\0", rbirif); + output->append_line(virtual_line); + sprintf(virtual_line, " (IM(FSAC(1,1))-IM(FSAC(2,2)))/IM(FSAC(1,1))=%15.7lE\n\0", rdichr); + output->append_line(virtual_line); + sprintf(virtual_line, " MULC\n\0"); + output->append_line(virtual_line); for (int i = 0; i < 4; i++) { - fprintf( - output, " %15.7lE%15.7lE%15.7lE%15.7lE\n", + sprintf( + virtual_line, " %15.7lE%15.7lE%15.7lE%15.7lE\n\0", cid->cmul[i][0], cid->cmul[i][1], cid->cmul[i][2], cid->cmul[i][3] ); + output->append_line(virtual_line); } - fprintf(output, " MULCLR\n"); + sprintf(virtual_line, " MULCLR\n\0"); + output->append_line(virtual_line); for (int i = 0; i < 4; i++) { - fprintf( - output, " %15.7lE%15.7lE%15.7lE%15.7lE\n", + sprintf( + virtual_line, " %15.7lE%15.7lE%15.7lE%15.7lE\n\0", cid->cmullr[i][0], cid->cmullr[i][1], cid->cmullr[i][2], cid->cmullr[i][3] ); + output->append_line(virtual_line); } if (iavm != 0) { mmulc(cid->c1ao->vintm, cid->cmullr, cid->cmul); @@ -1231,36 +1506,46 @@ int cluster_jxi488_cycle(int jxi488, ScattererConfiguration *sconf, GeometryConf for (int i = 0; i < 16; i++) { double value; value = real(cid->c1ao->vintm[i]); - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); value = imag(cid->c1ao->vintm[i]); - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); } for (int i = 0; i < 4; i++) { for (int j = 0; j < 4; j++) { double value = cid->cmul[i][j]; - tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + // tppoan.write(reinterpret_cast<char *>(&value), sizeof(double)); + vtppoanp->append_line(VirtualBinaryLine(value)); } } - fprintf(output, " CLUSTER (ENSEMBLE AVERAGE, MODE%2d)\n", iavm); + sprintf(virtual_line, " CLUSTER (ENSEMBLE AVERAGE, MODE%2d)\n\0", iavm); + output->append_line(virtual_line); if (inpol == 0) { - fprintf(output, " LIN\n"); + sprintf(virtual_line, " LIN\n\0"); + output->append_line(virtual_line); } else { // label 316 - fprintf(output, " CIRC\n"); + sprintf(virtual_line, " CIRC\n\0"); + output->append_line(virtual_line); } // label 318 - fprintf(output, " MULC\n"); + sprintf(virtual_line, " MULC\n\0"); + output->append_line(virtual_line); for (int i = 0; i < 4; i++) { - fprintf( - output, " %15.7lE%15.7lE%15.7lE%15.7lE\n", + sprintf( + virtual_line, " %15.7lE%15.7lE%15.7lE%15.7lE\n\0", cid->cmul[i][0], cid->cmul[i][1], cid->cmul[i][2], cid->cmul[i][3] ); + output->append_line(virtual_line); } - fprintf(output, " MULCLR\n"); + sprintf(virtual_line, " MULCLR\n\0"); + output->append_line(virtual_line); for (int i = 0; i < 4; i++) { - fprintf( - output, " %15.7lE%15.7lE%15.7lE%15.7lE\n", + sprintf( + virtual_line, " %15.7lE%15.7lE%15.7lE%15.7lE\n\0", cid->cmullr[i][0], cid->cmullr[i][1], cid->cmullr[i][2], cid->cmullr[i][3] ); + output->append_line(virtual_line); } } // label 420, continues jphs loop diff --git a/src/include/file_io.h b/src/include/file_io.h index e5c990083c4022c87e085c89cad2e8db213e9c42..df475bf7a0c7a2f1dac5f07cfb9d60c71ed363ba 100644 --- a/src/include/file_io.h +++ b/src/include/file_io.h @@ -5,6 +5,10 @@ #ifndef INCLUDE_FILE_IO_H_ #define INCLUDE_FILE_IO_H_ +#include <vector> + +class mixMPI; + /*! \class FileSchema * * \brief File content descriptor. @@ -156,4 +160,249 @@ class HDFFile { hid_t dapl_id=H5P_DEFAULT, hid_t dxpl_id=H5P_DEFAULT ); }; + +/*! \class VirtualAsciiFile + * + * \brief Virtual representation of an ASCII file. + */ +class VirtualAsciiFile { +protected: + //! \brief The number of lines. + // int32_t _num_lines; + //! \brief A vector of strings representing the file lines. + std::vector<std::string> *_file_lines; + +public: + // const int32_t &num_lines = _num_lines; + /*! \brief VirtualAsciiFile instance constructor. + * + * \param lines: `int32_t` Number of lines, if known in advance (optional, default is 0). + */ + VirtualAsciiFile(int32_t lines = 0); + + /*! \brief VirtualAsciiFile copy constructor. + * + * \param rhs: `const VirtualAsciiFile&` Reference to a VirtualAsciiFile instance. + */ + VirtualAsciiFile(const VirtualAsciiFile& rhs); + + /*! \brief VirtualAsciiFile instance constructor copying all contents off MPISend() calls from MPI process rr. + * + * \param mpidata: `mixMPI *` pointer to MPI data structure. + * \param rr: `int` rank of the MPI process sending the data. + */ + VirtualAsciiFile(const mixMPI *mpidata, int rr); + + /*! \brief VirtualAsciiFile instance destroyer. + */ + ~VirtualAsciiFile(); + + /*! \brief Append another VirtualAsciiFile at the end of the current instance. + * + * \param rhs: `const VirtualAsciiFile&` Reference to the VirtualAsciiFile to be appended. + */ + void append(VirtualAsciiFile& rhs); + + /*! \brief Append a line at the end of the file. + * + * \param line: `const string&` Reference to a string representing the line. + */ + void append_line(const std::string& line); + + /*! \brief Append the contents of the VirtualAsciiFile to a physical file on disk. + * + * \param file_name: `const string&` Name of the file to append contents to. + * \return result: `int` A result code (0 if successful). + */ + int append_to_disk(const std::string& file_name); + + /*! \brief Insert another VirtualAsciiFile at a given position. + * + * This function inserts a target VirtualAsciiFile in the current one at the given + * position. Optionally, a range of lines to be inserted can be specified, otherwise + * the full content of the target file is inserted. This function DOES NOT increase + * the size of the inner storage and it can only be used if the inner storage has + * already been adjusted to contain the insertion target. + * + * \param position: `int32_t` The position at which the other file is inserted in this one. + * \param rhs: `const VirtualAsciiFile&` The refence to the VirtualAsciiFile to be inserted. + * \param start: `int32_t` The first line to be inserted (optional, default is 0). + * \param end: `int32_t` The last line to be inserted (optional, default is 0 to read all). + * \param line: `const string&` Reference to a string representing the line. + * \return result: `int` A result code (0 if successful). + */ + int insert(int32_t position, VirtualAsciiFile& rhs, int32_t start = 0, int32_t end = 0); + + /*! \brief Get the number of lines in the current instance. + * + * \return size: `int32_t` The number of lines in the VirtualAsciiFile instance. + */ + int32_t number_of_lines() { return _file_lines->size(); } + + /*! \brief Write virtual file contents to a real file on disk. + * + * \param file_name: `const string&` Name of the file to append contents to. + * \return result: `int` A result code (0 if successful). + */ + int write_to_disk(const std::string& file_name); + + /*! \brief Send VirtualAsciiFile instance to MPI process 0 via MPISend() calls. + * + * \param mpidata: `mixMPI *` pointer to MPI data structure. + */ + void mpisend(const mixMPI *mpidata); +}; + + +/*! \class VirtualBinaryLine + * + * \brief Virtual representation of a binary file line. + */ +class VirtualBinaryLine { +// protected: +// //! \brief The pointer to the piece of data to be written, cast to char * +// char *_data_pointer; +// //! \brief the size of the data block. +// size_t _data_size; + +public: + //! \brief The pointer to the piece of data to be written, cast to char * + char *_data_pointer; + //! \brief the size of the data block. + size_t _data_size; + //! \brief Read only view of `_data_pointer`. + const char* data_pointer = _data_pointer; + //! \brief Read only view of `_data_size`. + const size_t & data_size = _data_size; + + /*! \brief VirtualBinaryLine instance constructor. + * + * \param mydata: `int, double, long, float, complex, or dcomplex`piece of data to put in the line. + */ + VirtualBinaryLine(int mydata); + VirtualBinaryLine(long mydata); + VirtualBinaryLine(float mydata); + VirtualBinaryLine(double mydata); + // VirtualBinaryLine(complex mydata); + VirtualBinaryLine(dcomplex mydata); + + /*! \brief VirtualBinaryLine copy constructor. + * + * \param rhs: `const VirtualBinaryLine&` Reference to a VirtualBinaryLine instance. + */ + VirtualBinaryLine(const VirtualBinaryLine& rhs); + + /*! \brief VirtualBinaryLine instance constructor copying all contents off MPISend() calls from MPI process rr. + * + * \param mpidata: `mixMPI *` pointer to MPI data structure. + * \param rr: `int` rank of the MPI process sending the data. + */ + VirtualBinaryLine(const mixMPI *mpidata, int rr); + + /*! \brief VirtualBinaryLine instance destroyer. + */ + ~VirtualBinaryLine(); + + /*! \brief Send VirtualBinaryLine instance to MPI process 0 via MPISend() calls. + * + * \param mpidata: `mixMPI *` pointer to MPI data structure. + */ + void mpisend(const mixMPI *mpidata); +}; + + +/*! \class VirtualBinaryFile + * + * \brief Virtual representation of a binary file. + */ +class VirtualBinaryFile { +protected: + //! \brief The number of lines. + // int32_t _num_lines; + // //! \brief A vector of strings representing the file lines. + // std::vector<VirtualBinaryLine> *_file_lines; + +public: + //! \brief A vector of strings representing the file lines. + std::vector<VirtualBinaryLine> *_file_lines; + // const int32_t &num_lines = _num_lines; + /*! \brief VirtualBinaryFile empty instance constructor. + * + */ + VirtualBinaryFile(); + + /*! \brief VirtualBinaryFile copy constructor. + * + * \param rhs: `const VirtualBinaryFile&` Reference to a VirtualBinaryFile instance. + */ + VirtualBinaryFile(const VirtualBinaryFile& rhs); + + /*! \brief VirtualBinaryFile instance constructor copying all contents off MPISend() calls from MPI process rr. + * + * \param mpidata: `mixMPI *` pointer to MPI data structure. + * \param rr: `int` rank of the MPI process sending the data. + */ + VirtualBinaryFile(const mixMPI *mpidata, int rr); + + /*! \brief VirtualBinaryFile instance destroyer. + */ + ~VirtualBinaryFile(); + + /*! \brief Append another VirtualBinaryFile at the end of the current instance. + * + * \param rhs: `const VirtualBinaryFile&` Reference to the VirtualBinaryFile to be appended. + */ + void append(VirtualBinaryFile& rhs); + + /*! \brief Append a line at the end of the file. + * + * \param line: `const string&` Reference to a string representing the line. + */ + void append_line(const VirtualBinaryLine& line); + + /*! \brief Append the contents of the VirtualBinaryFile to a physical file on disk. + * + * \param file_name: `const string&` Name of the file to append contents to. + * \return result: `int` A result code (0 if successful). + */ + int append_to_disk(const std::string& file_name); + + // /*! \brief Insert another VirtualBinaryFile at a given position. + // * + // * This function inserts a target VirtualBinaryFile in the current one at the given + // * position. Optionally, a range of lines to be inserted can be specified, otherwise + // * the full content of the target file is inserted. This function DOES NOT increase + // * the size of the inner storage and it can only be used if the inner storage has + // * already been adjusted to contain the insertion target. + // * + // * \param position: `int32_t` The position at which the other file is inserted in this one. + // * \param rhs: `const VirtualBinaryFile&` The refence to the VirtualBinaryFile to be inserted. + // * \param start: `int32_t` The first line to be inserted (optional, default is 0). + // * \param end: `int32_t` The last line to be inserted (optional, default is 0 to read all). + // * \param line: `const string&` Reference to a string representing the line. + // * \return result: `int` A result code (0 if successful). + // */ + // int insert(int32_t position, VirtualBinaryFile& rhs, int32_t start = 0, int32_t end = 0); + + /*! \brief Get the number of lines in the current instance. + * + * \return size: `int32_t` The number of lines in the VirtualBinaryFile instance. + */ + int32_t number_of_lines() { return _file_lines->size(); } + + /*! \brief Write virtual file contents to a real file on disk. + * + * \param file_name: `const string&` Name of the file to append contents to. + * \return result: `int` A result code (0 if successful). + */ + int write_to_disk(const std::string& file_name); + + /*! \brief Send VirtualBinaryFile instance to MPI process 0 via MPISend() calls. + * + * \param mpidata: `mixMPI *` pointer to MPI data structure. + */ + void mpisend(const mixMPI *mpidata); +}; #endif + + diff --git a/src/libnptm/file_io.cpp b/src/libnptm/file_io.cpp index c196d515e3923a5045373e28ce144c5a1e5de4c1..00afb19918befda10d592e2899ed2323a27a3811 100644 --- a/src/libnptm/file_io.cpp +++ b/src/libnptm/file_io.cpp @@ -3,7 +3,9 @@ * \brief Implementation of file I/O operations. */ #include <exception> +#include <fstream> #include <regex> +#include <cstring> #include <string> #include <hdf5.h> @@ -15,12 +17,21 @@ #include "../include/List.h" #endif +#ifndef INCLUDE_TYPES_H_ +#include "../include/types.h" +#endif + #ifndef INCLUDE_FILE_IO_H_ #include "../include/file_io.h" #endif +#ifdef USE_MPI +#include <mpi.h> +#endif + using namespace std; +/* >>> FileSchema class implementation <<< */ FileSchema::FileSchema(int num_rec, const std::string *rec_types, const std::string *rec_names) { num_records = num_rec; record_types = new string[num_rec]; @@ -48,7 +59,9 @@ string* FileSchema::get_record_types() { for (int i = 0; i < num_records; i++) rec_types[i] = record_types[i]; return rec_types; } +/* >>> End of FileSchema class implementation <<< */ +/* >>> HDFFile class implementation <<< */ HDFFile::HDFFile(const std::string& name, unsigned int flags, hid_t fcpl_id, hid_t fapl_id) { file_name = name; if (flags == H5F_ACC_EXCL || flags == H5F_ACC_TRUNC) @@ -222,3 +235,336 @@ herr_t HDFFile::write( } return status; } +/* >>> End of HDFFile class implementation <<< */ + +/* >>> VirtualAsciiFile class implementation <<< */ +VirtualAsciiFile::VirtualAsciiFile(int32_t lines) { + _file_lines = new vector<string>(); + for (int32_t li = 0; li < lines; li++) { + _file_lines->push_back(""); + } +} + +VirtualAsciiFile::VirtualAsciiFile(const VirtualAsciiFile& rhs) { + // _num_lines = rhs._num_lines; + _file_lines = new vector<string>(); + for (vector<string>::iterator it = rhs._file_lines->begin(); it != rhs._file_lines->end(); ++it) { + _file_lines->push_back(*it); + } +} + +#ifdef MPI_VERSION +VirtualAsciiFile::VirtualAsciiFile(const mixMPI *mpidata, int rr) { + // receive _num_lines from MPI process rr + int32_t num_lines; + MPI_Recv(&num_lines, 1, MPI_INT32_T, rr, 10, MPI_COMM_WORLD, MPI_STATUS_IGNORE); + int32_t mysize; + _file_lines = new vector<string>(); + // loop over data to receive + for (int32_t zi=0; zi<num_lines; zi++) { + // receive the size of the string to receive + MPI_Recv(&mysize, 1, MPI_INT32_T, rr, 10, MPI_COMM_WORLD, MPI_STATUS_IGNORE); + // allocate the buffer accordingly + char buffer[mysize+1]; + // receive the char buffer + MPI_Recv(buffer, mysize+1, MPI_CHAR, rr, 10, MPI_COMM_WORLD, MPI_STATUS_IGNORE); + // append to the vector + _file_lines->push_back(buffer); + } +} +#endif + +VirtualAsciiFile::~VirtualAsciiFile() { + // is it necessary to pop them out one by one? isn't there the dedicated method of std::vector to clean the vector? + // besides, shouldn't this be done anyway by the destructor of std:vector? + while (!_file_lines->size() > 0) { + _file_lines->pop_back(); + } + if (_file_lines != NULL) delete _file_lines; +} + +void VirtualAsciiFile::append(VirtualAsciiFile& rhs) { + // concatenate the virtualasciifile pointed by rhs to the current one + // can't we use the dedicated method insert of std::vector to do the appending, instead of an explicit loop? + for (vector<string>::iterator it = rhs._file_lines->begin(); it != rhs._file_lines->end(); ++it) { + _file_lines->push_back(*it); + } +} + +void VirtualAsciiFile::append_line(const string& line) { + // would it be worth reimplementing a sprintf-like method, so that we can give it all the arguments we would give to sprintf and get rid of the intermediate buffer completely? + // append a line of output to the virtualasciifile + _file_lines->push_back(line); +} + +int VirtualAsciiFile::append_to_disk(const std::string& file_name) { + // dump to disk the contents of the virtualasciifile, appending at the end of the given file_name + int result = 0; + fstream output_file; + output_file.open(file_name, ios::app); + if (output_file.is_open()) { + for (vector<string>::iterator it = _file_lines->begin(); it != _file_lines->end(); ++it) { + output_file << *it; + } + output_file.close(); + } else { + result = 1; + } + return result; +} + +int VirtualAsciiFile::insert(int32_t position, VirtualAsciiFile& rhs, int32_t start, int32_t end) { + int result = 0; + if (start == 0 && end == 0) { + end = rhs.number_of_lines(); + } + int32_t final_index = position + end - start; + if (final_index <= number_of_lines()) { + for (int32_t li = start; li < end; li++) { + // since here we are replacing the previous placeholder empty strings, make sure they are properly released when they are replaced (i.e. try it with a simple hello world example and pass it through valgrind) + _file_lines->at(position++) = rhs._file_lines->at(li); + } + } else { + // ERROR: target file is too long; + result = 1; + } + return result; +} + +int VirtualAsciiFile::write_to_disk(const std::string& file_name) { + // dump to disk the contents of the virtualasciifile, replacing the given file_name + int result = 0; + fstream output_file; + output_file.open(file_name, ios::out); + if (output_file.is_open()) { + for (vector<string>::iterator it = _file_lines->begin(); it != _file_lines->end(); ++it) { + output_file << *it; + } + output_file.close(); + } else { + result = 1; + } + return result; +} + +#ifdef MPI_VERSION +void VirtualAsciiFile::mpisend(const mixMPI *mpidata) { + // Send VirtualAsciiFile instance to MPI process 0 via MPISend() calls + // first send the size + int32_t num_lines = _file_lines->size(); + MPI_Send(&num_lines, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD); + // now loop over data to send + int32_t mysize; + for (vector<string>::iterator it = _file_lines->begin(); it != _file_lines->end(); ++it) { + // send the length of each string + mysize = (int32_t) it->size(); + MPI_Send(&mysize, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD); + // send the string itself + MPI_Send(it->c_str(), mysize+1, MPI_CHAR, 0, 10, MPI_COMM_WORLD); + } +} +#endif + +/* >>> End of VirtualAsciiFile class implementation <<< */ + +// /* >>> VirtualBinaryLine class implementation <<< */ +VirtualBinaryLine::VirtualBinaryLine(int mydata) { + _data_size = sizeof(mydata); + int *buffer = (int *) malloc(_data_size); + *buffer = mydata; + _data_pointer = reinterpret_cast<char *>(buffer); +} + +VirtualBinaryLine::VirtualBinaryLine(double mydata) { + _data_size = sizeof(mydata); + double *buffer = (double *) malloc(_data_size); + *buffer = mydata; + _data_pointer = reinterpret_cast<char *>(buffer); +} + +VirtualBinaryLine::VirtualBinaryLine(float mydata) { + _data_size = sizeof(mydata); + float *buffer = (float *) malloc(_data_size); + *buffer = mydata; + _data_pointer = reinterpret_cast<char *>(buffer); +} + +VirtualBinaryLine::VirtualBinaryLine(long mydata) { + _data_size = sizeof(mydata); + long *buffer = (long *) malloc(_data_size); + *buffer = mydata; + _data_pointer = reinterpret_cast<char *>(buffer); +} + +VirtualBinaryLine::VirtualBinaryLine(dcomplex mydata) { + _data_size = sizeof(mydata); + dcomplex *buffer = (dcomplex *) malloc(_data_size); + *buffer = mydata; + _data_pointer = reinterpret_cast<char *>(buffer); +} + +// VirtualBinaryLine::VirtualBinaryLine(complex mydata) { +// _data_size = sizeof(mydata); +// complex *buffer = (complex *) malloc(_data_size); +// *buffer = mydata; +// _data_pointer = reinterpret_cast<char *>(buffer); +// } + +VirtualBinaryLine::VirtualBinaryLine(const VirtualBinaryLine& rhs) { + _data_size = rhs._data_size; + _data_pointer = reinterpret_cast<char *>(malloc(rhs._data_size)); + memcpy(_data_pointer, rhs._data_pointer, _data_size); +} + +#ifdef MPI_VERSION +VirtualBinaryLine::VirtualBinaryLine(const mixMPI *mpidata, int rr) { + // receive mysize from MPI process rr + int32_t mysize; + MPI_Recv(&mysize, 1, MPI_INT32_T, rr, 10, MPI_COMM_WORLD, MPI_STATUS_IGNORE); + _data_size = mysize; + // allocate the buffer accordingly + _data_pointer = reinterpret_cast<char *>(malloc(mysize)); + // receive the char buffer + MPI_Recv(_data_pointer, mysize, MPI_CHAR, rr, 10, MPI_COMM_WORLD, MPI_STATUS_IGNORE); +} + +#endif + +VirtualBinaryLine::~VirtualBinaryLine() { + if (_data_pointer != NULL) { + free(_data_pointer); + _data_pointer = NULL; + } +} + +#ifdef MPI_VERSION +void VirtualBinaryLine::mpisend(const mixMPI *mpidata) { + // Send VirtualBinaryLine instance to MPI process 0 via MPISend() calls + // first send the size + int32_t mysize = _data_size; + MPI_Send(&mysize, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD); + // now send the data + MPI_Send(_data_pointer, mysize, MPI_CHAR, 0, 10, MPI_COMM_WORLD); +} +#endif +/* >>> End of VirtualBinaryLine class implementation <<< */ + + +/* >>> VirtualBinaryFile class implementation <<< */ +VirtualBinaryFile::VirtualBinaryFile() { + _file_lines = new vector<VirtualBinaryLine>(); +} + +VirtualBinaryFile::VirtualBinaryFile(const VirtualBinaryFile& rhs) { + // _num_lines = rhs._num_lines; + _file_lines = new vector<VirtualBinaryLine>(); + for (vector<VirtualBinaryLine>::iterator it = rhs._file_lines->begin(); it != rhs._file_lines->end(); ++it) { + _file_lines->push_back(VirtualBinaryLine(*it)); + } +} + +#ifdef MPI_VERSION +VirtualBinaryFile::VirtualBinaryFile(const mixMPI *mpidata, int rr) { + // receive _num_lines from MPI process rr + int32_t num_lines; + MPI_Recv(&num_lines, 1, MPI_INT32_T, rr, 10, MPI_COMM_WORLD, MPI_STATUS_IGNORE); + _file_lines = new vector<VirtualBinaryLine>(); + // loop over data to receive + for (int32_t zi=0; zi<num_lines; zi++) { + // receive the line of data + _file_lines->push_back(VirtualBinaryLine(mpidata, rr)); + } +} +#endif + +VirtualBinaryFile::~VirtualBinaryFile() { + // is it necessary to pop them out one by one? isn't there the dedicated method of std::vector to clean the vector? + // besides, shouldn't this be done anyway by the destructor of std:vector? + // for (vector<VirtualBinaryLine>::iterator it = _file_lines->begin(); it != _file_lines->end(); ++it) { + // delete it; + // } + while (!_file_lines->size() > 0) { + _file_lines->pop_back(); + } + if (_file_lines != NULL) delete _file_lines; +} + +void VirtualBinaryFile::append(VirtualBinaryFile& rhs) { + // concatenate the virtualasciifile pointed by rhs to the current one + // can't we use the dedicated method insert of std::vector to do the appending, instead of an explicit loop? + for (vector<VirtualBinaryLine>::iterator it = rhs._file_lines->begin(); it != rhs._file_lines->end(); ++it) { + _file_lines->push_back(VirtualBinaryLine(*it)); + } +} + +void VirtualBinaryFile::append_line(const VirtualBinaryLine& line) { + // would it be worth reimplementing a sprintf-like method, so that we can give it all the arguments we would give to sprintf and get rid of the intermediate buffer completely? + // append a line of output to the virtualasciifile + _file_lines->push_back(VirtualBinaryLine(line)); +} + +int VirtualBinaryFile::append_to_disk(const std::string& file_name) { + // dump to disk the contents of the virtualasciifile, appending at the end of the given file_name + int result = 0; + fstream output_file; + output_file.open(file_name, ios::app | ios::binary); + if (output_file.is_open()) { + for (vector<VirtualBinaryLine>::iterator it = _file_lines->begin(); it != _file_lines->end(); ++it) { + output_file.write(it->_data_pointer, it->_data_size); + } + output_file.close(); + } else { + result = 1; + } + return result; +} + +// int VirtualBinaryFile::insert(int32_t position, VirtualBinaryFile& rhs, int32_t start, int32_t end) { +// int result = 0; +// if (start == 0 && end == 0) { +// end = rhs.number_of_lines(); +// } +// int32_t final_index = position + end - start; +// if (final_index <= number_of_lines()) { +// for (int32_t li = start; li < end; li++) { +// // since here we are replacing the previous placeholder empty strings, make sure they are properly released when they are replaced (i.e. try it with a simple hello world example and pass it through valgrind) +// VirtualBinaryLine templine = VirtualBinaryLine(rhs._file_lines->at(li)); +// _file_lines->at(position++) = templine; +// } +// } else { +// // ERROR: target file is too long; +// result = 1; +// } +// return result; +// } + +int VirtualBinaryFile::write_to_disk(const std::string& file_name) { + // dump to disk the contents of the virtualasciifile, replacing the given file_name + int result = 0; + fstream output_file; + output_file.open(file_name, ios::out | ios::binary); + if (output_file.is_open()) { + for (vector<VirtualBinaryLine>::iterator it = _file_lines->begin(); it != _file_lines->end(); ++it) { + output_file.write(it->_data_pointer, it->_data_size); + } + output_file.close(); + } else { + result = 1; + } + return result; +} + +#ifdef MPI_VERSION +void VirtualBinaryFile::mpisend(const mixMPI *mpidata) { + // Send VirtualBinaryFile instance to MPI process 0 via MPISend() calls + // first send the size + int32_t num_lines = _file_lines->size(); + MPI_Send(&num_lines, 1, MPI_INT32_T, 0, 10, MPI_COMM_WORLD); + // now loop over data to send + for (vector<VirtualBinaryLine>::iterator it = _file_lines->begin(); it != _file_lines->end(); ++it) { + it->mpisend(mpidata); + } +} +#endif + +/* >>> End of VirtualBinaryFile class implementation <<< */ diff --git a/src/make.inc b/src/make.inc index 6f65edda86bf6399edff50f2172a70007a2d7e5d..57499cbead9efdd7acd6ad3658982013e8b69df9 100644 --- a/src/make.inc +++ b/src/make.inc @@ -124,7 +124,7 @@ endif # CXXFLAGS defines the default compilation options for the C++ compiler ifndef CXXFLAGS - override CXXFLAGS=-O3 -ggdb -pg -coverage -I$(HDF5_INCLUDE) $(MPI_CXXFLAGS) $(NVTX_CXXFLAGS) + override CXXFLAGS=-O3 -ggdb -pg -coverage -Wno-format-contains-nul -I$(HDF5_INCLUDE) $(MPI_CXXFLAGS) $(NVTX_CXXFLAGS) ifdef USE_OPENMP override CXXFLAGS+= -fopenmp # closes USE_OPENMP diff --git a/src/scripts/pycompare.py b/src/scripts/pycompare.py index 193860f0b95245e47976bc1be916374d75e3837d..a5129deccf54c1071bd2c4d6f002defe18a9ffa4 100755 --- a/src/scripts/pycompare.py +++ b/src/scripts/pycompare.py @@ -115,6 +115,8 @@ def compare_files(config): print("INFO: using line-wise mode") print("INFO: counting result lines...") while (f_lines[0] != ''): + if (c_lines[0].startswith("INSERTION:")): + c_lines = [c_file.readline()] if (c_lines[0] != ''): line_count += 1 else: @@ -166,6 +168,8 @@ def compare_files(config): else: f_lines = [fortran_file.readline()] c_lines = [c_file.readline()] + if (c_lines[0].startswith("INSERTION:")): + c_lines = [c_file.readline()] num_read_lines += 1 # Start here the comparison loop if (len(f_lines) == len(c_lines)): @@ -429,7 +433,7 @@ def parse_arguments(): 'fortran_file_name': '', 'c_file_name': '', 'full_log': False, - 'linewise': False, + 'linewise': True, 'log_html': False, 'html_output': 'pycompare.html', 'warning_threshold': 0.005, @@ -475,7 +479,7 @@ def print_help(): print("--full Print all lines to log file (default prints only mismatches).") print("--help Print this help and exit.") print("--html[=OPT_OUTPUT_NAME] Enable logging to HTML file (default logs to \"pycompare.html\").") - print("--linewise Load only one line at a time. Useful to compare big files (false by default).") + print("--linewise Load only one line at a time. Useful to compare big files (True by default).") print("--quick Stop on first mismatch (default is to perform a full check).") print("--warn Set a fractional threshold for numeric warning (default = 0.005).") print(" ")