00001
00017
00018 #include "femaxx_post3d.hh"
00019 #include "femaxx_post3dparam.hh"
00020
00021
00022 int main(int argc, char** argv)
00023 {
00024
00025
00026
00027 #ifdef HAVE_RLOG
00028
00029 rlog::StdioNode stdLog;
00030
00034 #ifdef USE_RLOG_DEBUG_CHANNEL
00035 stdLog.subscribeTo(rlog::GetGlobalChannel("debug"));
00036 #endif //USE_RLOG_DEBUG_CHANNEL
00037 #ifdef USE_RLOG_ERROR_CHANNEL
00038 stdLog.subscribeTo(rlog::GetGlobalChannel("error"));
00039 #endif //USE_RLOG_ERROR_CHANNEL
00040 #ifdef USE_RLOG_INFO_CHANNEL
00041 stdLog.subscribeTo(rlog::GetGlobalChannel("info"));
00042 #endif //USE_RLOG_INFO_CHANNEL
00043 #ifdef USE_RLOG_WARNING_CHANNEL
00044 stdLog.subscribeTo(rlog::GetGlobalChannel("warning"));
00045 #endif //USE_RLOG_WARNING_CHANNEL
00046 #endif // HAVE_RLOG
00047
00048
00049 #ifdef HAVE_RLOG
00050 rInfo(" ");
00051 rInfo("this is femaxx_post3d");
00052 rInfo(" ");
00053 rInfo("the program converts femaxx computational results into VTK files and ANSYS readable files");
00054 rInfo("2006-2007, all rights reserved, benedikt oswald");
00055 rInfo("http://maxwell.psi.ch/amaswiki/index.php/User:BenediktOswald");
00056 rInfo(" ");
00057 rInfo("history: 2006 jan 27, benedikt oswald, creation");
00058 rInfo(" 2006 feb 10, benedikt oswald, reorganized into main and several functions");
00059 rInfo(" 2007 jun 04, benedikt oswald, added logging functionality and command line parameter processing ");
00060 rInfo("capabilities : read results femaxx eigenmodal calculations results from hdf5 file");
00061 rInfo(" and converts them into VTK files and files that can be read by ANSYS");
00062 rInfo(" for subsequent further analsys.");
00063 #else
00064 std::cout << "this is femaxx_post3d" << std::endl;
00065 std::cout << "the program converts femaxx computational results into VTK files and ANSYS readable files" << std::endl;
00066 std::cout << "2006-2007, all rights reserved, benedikt oswald" << std::endl;
00067 #endif // HAVE_RLOG
00068
00069
00070
00071
00072
00073
00075 FemaxxPost3dparam *commandLineParameter;
00076 commandLineParameter = FemaxxPost3dparam::getInstance();
00077
00078
00079 try
00080 {
00081 std::vector<unsigned int> startValueDef(1,0);
00083 boost::program_options::options_description
00084 mandatory_options("mandatory program options");
00085 mandatory_options.add_options()
00086 ("help,h", "produce this help")
00087 ("input-file", boost::program_options::value<std::string>(),
00088 "this is the HDF5 file generated by femaxx, containing the eigenmodal solution.")
00089 ("output-file", boost::program_options::value<std::string>(),
00090 "this is the base file name for the VTK and ANSYS readable output files.");
00091
00093 boost::program_options::options_description
00094 tag_mapping_options("optional program Options");
00095
00096
00097
00098
00099
00100
00101
00102
00103
00104
00105
00106
00107
00108
00109
00110
00111
00113 boost::program_options::options_description options;
00114 options.add(mandatory_options).add(tag_mapping_options);
00115
00117 boost::program_options::variables_map varMap;
00118 boost::program_options::store
00119 (boost::program_options::parse_command_line(argc, argv, options),
00120 varMap);
00121 boost::program_options::notify(varMap);
00122
00124 if (varMap.count("help"))
00125 {
00127 cout << options << "\n";
00128 return(OKCODE);
00129 }
00130 else
00131 {
00133 if (varMap.count("input-file"))
00134 commandLineParameter->inputFile(varMap["input-file"].as<std::string>());
00135 if (varMap.count("output-file"))
00136 commandLineParameter->outputFile(varMap["output-file"].as<std::string>());
00137
00138
00139
00140
00141
00142
00143
00144 }
00145 }
00147 catch(exception& error)
00148 {
00149 #ifdef HAVE_RLOG
00150 rError("Unknown command line option.");
00151 rError("Use '--help' for detail informations.");
00152 rError("Error code from boost programm_options: %s",error.what());
00153 rError("Exit!");
00154 #else
00155 std::cerr << "Unknown command line option." << std::endl
00156 << "Use '--help' for detail informations." << std::endl
00157 << "Error code from boost programm_options: %s"
00158 << error.what() << std::endl
00159 << "Exit!" << std::endl;
00160 #endif // HAVE_RLOG
00161 return(ERRORCODE);
00162 }
00163
00164
00165
00166
00167
00168 #ifdef HAVE_RLOG
00169 rInfo("============================================================================");
00170 rInfo("femaxx eigenmodal calculation results in HDF5 file to VTK legacy file format");
00171 rInfo("============================================================================");
00172 rInfo(" ");
00173 #else
00174 std::cout
00175 << "============================================================================"
00176 << std::endl
00177 << "femaxx eigenmodal calculation results in HDF5 file to VTK legacy file format"
00178 << std::endl
00179 << "============================================================================"
00180 << std::endl;
00181 #endif // HAVE_RLOG
00182
00183
00185 commandLineParameter->showAll();
00187 if (commandLineParameter->checkAll() == (ERRORCODE) )
00188 {
00189 #ifdef HAVE_RLOG
00190 rError("Some mandatory command line parameters are not specified.");
00191 rError("Exit!");
00192 #else
00193 std::cerr << "Some mandatory command line parameters are not specified."
00194 << std::endl << "Error! " << std::endl;
00195 #endif // HAVE_RLOG
00196 return(ERRORCODE);
00197 }
00198
00199
00202
00203 int m;
00204 int res;
00205 int npoint;
00206 int ntet;
00207 int nsample;
00208
00209 id_t** tetvertid=NULL;
00210 double** coord3d=NULL;
00211
00212
00213 std::string eigenfile(commandLineParameter->inputFile());
00214 #ifdef HAVE_RLOG
00215 rInfo("reading eigenmodal solution from file : %s", commandLineParameter->inputFile().c_str());
00216 #endif // HAVE_RLOG
00217
00218
00219
00220
00221
00224
00225
00226
00227
00229 #ifdef HAVE_RLOG
00230 rInfo("..reading vertex and mesh information from hdf5 file");
00231 #endif // HAVE_RLOG
00232
00233
00234 res=h5_read_point(eigenfile, H5_GROUP_FEMAX_EIGENMESH, npoint, coord3d);
00235 res=h5_read_tetvertid(eigenfile, H5_GROUP_FEMAX_EIGENMESH, ntet, tetvertid);
00236
00237 #ifdef HAVE_RLOG
00238 rInfo("read eigenmodal solution data from hdf5 file");
00239 rInfo("eigenmodal solution: mesh has %d tetrahedra", ntet);
00240 rInfo("eigenmodal solution: mesh has %d vertices", npoint);
00241 #endif // HAVE_RLOG
00242
00244 #ifdef HAVE_RLOG
00245 rInfo("...writing mesh in legacy VTK file format");
00246 #endif // HAVE_RLOG
00247
00248 std::string vtklffmesh(commandLineParameter->outputFile());
00249 vtklffmesh.append(".mesh");
00250 res=vtk_write_mesh(vtklffmesh, npoint, coord3d, ntet, tetvertid);
00251
00252 #ifdef HAVE_RLOG
00253 rInfo("wrote mesh in legacy VTK file format");
00254 #endif // HAVE_RLOG
00255
00256
00257
00258 #ifdef HAVE_RLOG
00259 rInfo("...reading eigenvalues, eigenfrequencies and quality factors from hdf5 file");
00260 #endif // HAVE_RLOG
00261
00262 double* lambda=NULL;
00263 int nmode=h5_read_eigenvalues(eigenfile,H5_GROUP_FEMAX_EIGENMODES,lambda);
00264
00266 vector<vector<double> > eigenmodalsolution(nmode);
00267
00268 unsigned int mode=0;
00269
00270 for(vector<vector<double> >::iterator it=eigenmodalsolution.begin(); it != eigenmodalsolution.end(); it++)
00271 {
00272 (*it).push_back(lambda[mode*3]);
00273 (*it).push_back(lambda[mode*3+1]);
00274 (*it).push_back(lambda[mode*3+2]);
00275 mode++;
00276 }
00277
00278
00279 #ifdef HAVE_RLOG
00280 rInfo("computed eigenmodal quantity:");
00281 rInfo("eigenvalue [-] :: resonance frequency [Hz] :: quality figure [-]");
00282 rInfo("---------------------------------------------------------------------------------------");
00283
00284 mode=0;
00285
00286 for(vector<vector<double> >::iterator it=eigenmodalsolution.begin(); it != eigenmodalsolution.end(); it++)
00287 {
00288 rInfo("mode [%d] : %12.9e \t %12.9e \t %12.9e", mode, (*it)[0], (*it)[1], (*it)[2] );
00289 mode++;
00290 }
00291
00292 rInfo("---------------------------------------------------------------------------------------");
00293
00294 #endif // HAVE_RLOG
00295
00296
00297
00298 #ifdef HAVE_RLOG
00299 rInfo("read eigenvalues, eigenfrequencies and quality factors from hdf5 file");
00300 #endif // HAVE_RLOG
00301
00302
00303
00304 #ifdef HAVE_RLOG
00305 rInfo("eigenmodal solution contains %d modes", nmode);
00306 #endif // HAVE_RLOG
00307
00308 double*** electricfields=new double**[nmode];
00309 double*** magneticfields=new double**[nmode];
00310
00311
00312 #ifdef HAVE_RLOG
00313 rInfo("eigenmodal solution contains %d modes", nmode);
00314 #endif // HAVE_RLOG
00315
00316
00317 for(m=0;m<nmode;m++)
00318 {
00319
00320
00321 #ifdef HAVE_RLOG
00322 rInfo(" ");
00323 rInfo("...reading eigenmodal solution no. %d [electric and magnetic field]", m);
00324 #endif // HAVE_RLOG
00325
00326 res=h5_read_sampledefield(eigenfile, H5_GROUP_FEMAX_EIGENMODES, nsample, electricfields[m], m);
00327 res=h5_read_sampledhfield(eigenfile, H5_GROUP_FEMAX_EIGENMODES, nsample, magneticfields[m], m);
00328
00329 #ifdef HAVE_RLOG
00330 rInfo("read eigenmodal solution [%d]", m);
00331 #endif // HAVE_RLOG
00332
00333 }
00334
00335 #ifdef HAVE_RLOG
00336 rInfo("read %d eigenmodal solution [electric and magnetic field]", nmode);
00337 rInfo(" ");
00338 #endif // HAVE_RLOG
00339
00341 double** sloc3d=NULL;
00342 int nsloc=0;
00343
00344 res=h5_read_sloc(eigenfile, H5_GROUP_FEMAX_EIGENMODES, nsloc, sloc3d);
00345
00346
00348 #ifdef HAVE_RLOG
00349 rInfo("...writing eigenmodal field solution in legacy VTK file format");
00350 #endif // HAVE_RLOG
00351
00352 std::string vtklffehfield(commandLineParameter->outputFile());
00353 vtklffehfield.append(".field");
00354
00355 res=vtk_write_field(vtklffehfield,
00356 npoint,
00357 coord3d,
00358 ntet,
00359 tetvertid,
00360 electricfields,
00361 magneticfields,
00362 nmode);
00363
00364 #ifdef HAVE_RLOG
00365 rInfo("wrote eigenmodal field solution in legacy VTK file format");
00366 #endif // HAVE_RLOG
00367
00368
00369
00370 #ifdef HAVE_RLOG
00371 rInfo("...writing eigenmodal field solution in ANSYS readable format");
00372 #endif // HAVE_RLOG
00373
00374 std::string scffield(eigenfile);
00375 res=scf_write_efield(scffield, nsample, sloc3d, electricfields, nmode);
00376 res=scf_write_hfield(scffield, nsample, sloc3d, magneticfields, nmode);
00377
00378 #ifdef HAVE_RLOG
00379 rInfo("wrote eigenmodal field solution in ANSYS readable format");
00380 #endif // HAVE_RLOG
00381
00382
00387
00388 return(OKCODE);
00389 }
00390
00391
00392
00394 int h5_read_point(std::string filename, std::string groupname, int &npoint, double** &coord3d)
00395 {
00396 int i,j;
00399 int rank;
00400 int npoint_=0;
00401 hsize_t dimsf[2];
00402 hid_t plist_id, file_id, group_id, dataset_id, dataspace_id, filespace_id, memspace_id;
00403 herr_t status;
00404
00405 plist_id = H5Pcreate(H5P_FILE_ACCESS);
00406 file_id = H5Fopen(filename.c_str(), H5F_ACC_RDONLY, plist_id);
00407 assert(file_id >= 0);
00408 H5Pclose(plist_id);
00409
00410 group_id = H5Gopen(file_id, groupname.c_str());
00411
00413 dataset_id = H5Dopen(group_id,H5_NAME_FEMAX_POINTS.c_str());
00414 dataspace_id = H5Dget_space(dataset_id);
00415 rank = H5Sget_simple_extent_ndims(dataspace_id);
00416 status = H5Sget_simple_extent_dims(dataspace_id, dimsf,NULL);
00417 filespace_id = H5Screate_simple(2, dimsf, NULL);
00418 npoint=dimsf[0];
00419 npoint_=npoint;
00420
00422 hsize_t offset[2];
00423 hsize_t count[2];
00424
00425 offset[0] = 0;
00426 offset[1] = 0;
00427 count[0] = npoint_;
00428 count[1] = NREALCOORD3D;
00429
00431 status = H5Sselect_hyperslab(dataspace_id, H5S_SELECT_SET, offset, NULL, count, NULL);
00432
00434 hsize_t dimsm[2];
00435
00436 dimsm[0]=dimsf[0];
00437 dimsm[1]=dimsf[1];
00438 memspace_id=H5Screate_simple(rank,dimsm,NULL);
00439
00441 hsize_t offsetm[2];
00442 hsize_t countm[2];
00443
00444 offsetm[0]=0;
00445 offsetm[1]=0;
00446
00447 countm[0]=dimsf[0];
00448 countm[1]=dimsf[1];
00449
00450 status=H5Sselect_hyperslab(memspace_id,H5S_SELECT_SET, offsetm, NULL, countm, NULL);
00451
00453 double* coord3dlin=new double[npoint_*NREALCOORD3D];
00454
00455
00457 status=H5Dread(dataset_id,H5T_NATIVE_DOUBLE,memspace_id,dataspace_id,H5P_DEFAULT, coord3dlin);
00458
00459 coord3d=new double*[npoint_];
00460 for(i=0;i<npoint_;i++)
00461 {
00462 coord3d[i]=new double[NREALCOORD3D];
00463 for(j=0;j<NREALCOORD3D;j++)
00464 {
00465 coord3d[i][j]=coord3dlin[i*NREALCOORD3D+j];
00466 }
00467 }
00468
00469 #ifdef __DEBUG__VERBOSE__
00470 cout << "[DEBUGOUT] : " << __FILE__ << " : " << __LINE__ << endl;
00471
00472 for(i=0;i<npoint_;i++)
00473 {
00474 for(j=0;j<NREALCOORD3D;j++)
00475 {
00476 cout << coord3d[i][j] << TAB_STRING.c_str();
00477 }
00478 cout << endl;
00479 }
00480
00481
00482 #endif
00485 H5Dclose(dataset_id);
00486 H5Sclose(dataspace_id);
00487 H5Sclose(memspace_id);
00488 H5Gclose(group_id);
00489 H5Fclose(file_id);
00490
00492 delete []coord3dlin;
00493
00494
00495 return(OKCODE);
00496 }
00497
00499 int h5_read_tetvertid(std::string filename, std::string groupname, int &ntet, id_t** &tetvertid)
00500 {
00501 int i,j;
00504 int rank;
00505 int ncol=0;
00506 int ntet_=0;
00507
00508 hsize_t dimsf[2];
00509 hsize_t offset[2];
00510 hsize_t count[2];
00511 hsize_t dimsm[2];
00512 hsize_t offsetm[2];
00513 hsize_t countm[2];
00514
00515 hid_t plist_id, file_id, group_id, dataset_id, dataspace_id, filespace_id, memspace_id;
00516 herr_t status;
00517
00518 plist_id = H5Pcreate(H5P_FILE_ACCESS);
00519 file_id = H5Fopen(filename.c_str(), H5F_ACC_RDONLY, plist_id);
00520 assert(file_id >= 0);
00521 H5Pclose(plist_id);
00522
00524 group_id = H5Gopen(file_id, groupname.c_str());
00525 dataset_id = H5Dopen(group_id,H5_NAME_FEMAX_TETVID.c_str());
00526 dataspace_id = H5Dget_space(dataset_id);
00527 rank = H5Sget_simple_extent_ndims(dataspace_id);
00528 status = H5Sget_simple_extent_dims(dataspace_id, dimsf,NULL);
00529 filespace_id = H5Screate_simple(2, dimsf, NULL);
00530 ntet=dimsf[0];
00531 ntet_=ntet;
00532 ncol=dimsf[1];
00533
00535 offset[0] = 0;
00536 offset[1] = 0;
00537 count[0] = ntet_;
00538 count[1] = ncol;
00539
00541 status = H5Sselect_hyperslab(dataspace_id, H5S_SELECT_SET, offset, NULL, count, NULL);
00542
00544 dimsm[0]=dimsf[0];
00545 dimsm[1]=dimsf[1];
00546 memspace_id=H5Screate_simple(rank,dimsm,NULL);
00547
00549 offsetm[0]=0;
00550 offsetm[1]=0;
00551
00552 countm[0]=dimsf[0];
00553 countm[1]=dimsf[1];
00554
00555 status=H5Sselect_hyperslab(memspace_id,H5S_SELECT_SET, offsetm, NULL, countm, NULL);
00556
00558 id_t* tetvertexidlin=new id_t[ntet_*ncol];
00559
00560
00562 status=H5Dread(dataset_id,H5T_NATIVE_UINT,memspace_id,dataspace_id,H5P_DEFAULT, tetvertexidlin);
00563
00565 tetvertid=new id_t*[ntet_];
00566 for(i=0;i<ntet_;i++)
00567 {
00568 tetvertid[i]=new id_t[ncol];
00569
00570 for(j=0;j<ncol;j++)
00571 {
00572 tetvertid[i][j]=tetvertexidlin[i*ncol+j];
00573 }
00574 }
00575
00576 #ifdef __DEBUG__VERBOSE__
00577 cout << "[DEBUGOUT] : " << __FILE__ << " : " << __LINE__ << endl;
00578 for(i=0;i<ntet_;i++)
00579 {
00580 for(j=0;j<NCORNERTET;j++)
00581 {
00582 cout << tetvertid[i][j] << " : ";
00583 }
00584 cout << endl;
00585 }
00586 #endif
00592 delete []tetvertexidlin;
00593
00595 H5Dclose(dataset_id);
00596 H5Sclose(dataspace_id);
00597 H5Sclose(memspace_id);
00598 H5Gclose(group_id);
00599 H5Fclose(file_id);
00600
00601
00602 return(OKCODE);
00603 }
00604
00605
00607 int h5_read_sloc(std::string filename, std::string groupname, int &npoint, double** &coord3d)
00608 {
00609 int i,j;
00612 int rank;
00613 int npoint_=0;
00614 hsize_t dimsf[2];
00615 hid_t plist_id, file_id, group_id, dataset_id, dataspace_id, filespace_id, memspace_id;
00616 herr_t status;
00617
00618 plist_id = H5Pcreate(H5P_FILE_ACCESS);
00619 file_id = H5Fopen(filename.c_str(), H5F_ACC_RDONLY, plist_id);
00620 assert(file_id >= 0);
00621 H5Pclose(plist_id);
00622
00623 group_id = H5Gopen(file_id, H5_GROUP_FEMAX_EIGENMODES.c_str());
00624
00626 dataset_id = H5Dopen(group_id,H5_NAME_FEMAX_SLOCBASENAME.c_str());
00627 dataspace_id = H5Dget_space(dataset_id);
00628 rank = H5Sget_simple_extent_ndims(dataspace_id);
00629 status = H5Sget_simple_extent_dims(dataspace_id, dimsf,NULL);
00630 filespace_id = H5Screate_simple(2, dimsf, NULL);
00631 npoint=dimsf[0];
00632 npoint_=npoint;
00633
00635 hsize_t offset[2];
00636 hsize_t count[2];
00637
00638 offset[0] = 0;
00639 offset[1] = 0;
00640 count[0] = npoint_;
00641 count[1] = NREALCOORD3D;
00642
00644 status = H5Sselect_hyperslab(dataspace_id, H5S_SELECT_SET, offset, NULL, count, NULL);
00645
00647 hsize_t dimsm[2];
00648
00649 dimsm[0]=dimsf[0];
00650 dimsm[1]=dimsf[1];
00651 memspace_id=H5Screate_simple(rank,dimsm,NULL);
00652
00654 hsize_t offsetm[2];
00655 hsize_t countm[2];
00656
00657 offsetm[0]=0;
00658 offsetm[1]=0;
00659
00660 countm[0]=dimsf[0];
00661 countm[1]=dimsf[1];
00662
00663 status=H5Sselect_hyperslab(memspace_id,H5S_SELECT_SET, offsetm, NULL, countm, NULL);
00664
00666 double* coord3dlin=new double[npoint_*NREALCOORD3D];
00667
00668
00670 status=H5Dread(dataset_id,H5T_NATIVE_DOUBLE,memspace_id,dataspace_id,H5P_DEFAULT, coord3dlin);
00671
00672 coord3d=new double*[npoint_];
00673 for(i=0;i<npoint_;i++)
00674 {
00675 coord3d[i]=new double[NREALCOORD3D];
00676 for(j=0;j<NREALCOORD3D;j++)
00677 {
00678 coord3d[i][j]=coord3dlin[i*NREALCOORD3D+j];
00679 }
00680 }
00681
00682 #ifdef __DEBUG__VERBOSE__
00683 cout << "[DEBUGOUT] : " << __FILE__ << " : " << __LINE__ << endl;
00684
00685 for(i=0;i<npoint_;i++)
00686 {
00687 for(j=0;j<NREALCOORD3D;j++)
00688 {
00689 cout << coord3d[i][j] << TAB_STRING.c_str();
00690 }
00691 cout << endl;
00692 }
00693
00694
00695 #endif
00698 H5Dclose(dataset_id);
00699 H5Sclose(dataspace_id);
00700 H5Sclose(memspace_id);
00701 H5Gclose(group_id);
00702 H5Fclose(file_id);
00703
00705 delete []coord3dlin;
00706
00707
00708 return(OKCODE);
00709
00710
00711 }
00712
00713
00715 int h5_read_sampledefield(std::string filename, std::string groupname, int &nsample, double** &efield, int mthmode)
00716 {
00717 int i,j;
00720 int rank;
00721 int ncol=0;
00722 int nsample_=0;
00723
00724 hsize_t dimsf[2];
00725 hsize_t offset[2];
00726 hsize_t count[2];
00727 hsize_t dimsm[2];
00728 hsize_t offsetm[2];
00729 hsize_t countm[2];
00730
00731 hid_t plist_id, file_id, group_id, dataset_id, dataspace_id, filespace_id, memspace_id;
00732 herr_t status;
00733
00735 plist_id = H5Pcreate(H5P_FILE_ACCESS);
00736 file_id = H5Fopen(filename.c_str(), H5F_ACC_RDONLY, plist_id);
00737 assert(file_id >= 0);
00738 H5Pclose(plist_id);
00739 group_id = H5Gopen(file_id, H5_GROUP_FEMAX_EIGENMODES.c_str());
00740
00741
00742 std::string mthelectricfieldname(H5_NAME_FEMAX_EIGENFIELD);
00743 std::ostringstream out;
00744 out << mthmode;
00745
00746 mthelectricfieldname.append(out.str());
00747
00748 dataset_id = H5Dopen(group_id,mthelectricfieldname.c_str());
00749
00750 dataspace_id = H5Dget_space(dataset_id);
00751 rank = H5Sget_simple_extent_ndims(dataspace_id);
00752 status = H5Sget_simple_extent_dims(dataspace_id, dimsf,NULL);
00753 filespace_id = H5Screate_simple(2, dimsf, NULL);
00754 nsample=dimsf[0];
00755 nsample_=nsample;
00756 ncol=dimsf[1];
00757
00759 offset[0] = 0;
00760 offset[1] = 0;
00761 count[0] = nsample_;
00762 count[1] = ncol;
00763
00765 status = H5Sselect_hyperslab(dataspace_id, H5S_SELECT_SET, offset, NULL, count, NULL);
00766
00768 dimsm[0]=dimsf[0];
00769 dimsm[1]=dimsf[1];
00770 memspace_id=H5Screate_simple(rank,dimsm,NULL);
00771
00773 offsetm[0]=0;
00774 offsetm[1]=0;
00775
00776 countm[0]=dimsf[0];
00777 countm[1]=dimsf[1];
00778
00779 status=H5Sselect_hyperslab(memspace_id,H5S_SELECT_SET, offsetm, NULL, countm, NULL);
00780
00782 double* ev=new double[nsample_*ncol];
00783 for(i=0;i<nsample_*ncol;i++){ ev[i]=0.0; }
00784
00786 status=H5Dread(dataset_id,H5T_NATIVE_DOUBLE,memspace_id,dataspace_id,H5P_DEFAULT, ev);
00787
00788 #ifdef __DEBUG__VERBOSE__
00789 cout << "[DEBUGOUT] : " << __FILE__ << " : " << __LINE__ << endl;
00790 #endif
00793 efield=new double*[nsample_];
00794
00795 for(i=0;i<nsample_;i++)
00796 {
00797 efield[i]=new double[ncol];
00798
00799 for(j=0;j<ncol;j++)
00800 {
00801 efield[i][j]=ev[i*ncol+j];
00802 }
00803 }
00804
00806 delete []ev;
00807
00809 H5Dclose(dataset_id);
00810 H5Sclose(dataspace_id);
00811 H5Sclose(memspace_id);
00812 H5Gclose(group_id);
00813 H5Fclose(file_id);
00814
00815
00816 return(OKCODE);
00817 }
00818
00819
00821 int h5_read_sampledhfield(std::string filename, std::string groupname, int &nsample, double** &hfield, int mthmode)
00822 {
00823 int i,j;
00826 int rank;
00827 int ncol=0;
00828 int nsample_=0;
00829
00830 hsize_t dimsf[2];
00831 hsize_t offset[2];
00832 hsize_t count[2];
00833 hsize_t dimsm[2];
00834 hsize_t offsetm[2];
00835 hsize_t countm[2];
00836
00837 hid_t plist_id, file_id, group_id, dataset_id, dataspace_id, filespace_id, memspace_id;
00838 herr_t status;
00839
00841 plist_id = H5Pcreate(H5P_FILE_ACCESS);
00842 file_id = H5Fopen(filename.c_str(), H5F_ACC_RDONLY, plist_id);
00843 assert(file_id >= 0);
00844 H5Pclose(plist_id);
00845 group_id = H5Gopen(file_id, H5_GROUP_FEMAX_EIGENMODES.c_str());
00846
00847
00848 std::string mthemagneticfieldname(H5_NAME_FEMAX_EIGENCURL);
00849 std::ostringstream out;
00850 out << mthmode;
00851
00852 mthemagneticfieldname.append(out.str());
00853
00854 dataset_id = H5Dopen(group_id,mthemagneticfieldname.c_str());
00855
00856 dataspace_id = H5Dget_space(dataset_id);
00857 rank = H5Sget_simple_extent_ndims(dataspace_id);
00858 status = H5Sget_simple_extent_dims(dataspace_id, dimsf,NULL);
00859 filespace_id = H5Screate_simple(2, dimsf, NULL);
00860 nsample=dimsf[0];
00861 nsample_=nsample;
00862 ncol=dimsf[1];
00863
00865 offset[0] = 0;
00866 offset[1] = 0;
00867 count[0] = nsample_;
00868 count[1] = ncol;
00869
00871 status = H5Sselect_hyperslab(dataspace_id, H5S_SELECT_SET, offset, NULL, count, NULL);
00872
00874 dimsm[0]=dimsf[0];
00875 dimsm[1]=dimsf[1];
00876 memspace_id=H5Screate_simple(rank,dimsm,NULL);
00877
00879 offsetm[0]=0;
00880 offsetm[1]=0;
00881
00882 countm[0]=dimsf[0];
00883 countm[1]=dimsf[1];
00884
00885 status=H5Sselect_hyperslab(memspace_id,H5S_SELECT_SET, offsetm, NULL, countm, NULL);
00886
00888 double* hv=new double[nsample_*ncol];
00889 for(i=0;i<nsample_*ncol;i++){ hv[i]=0.0; }
00890
00892 status=H5Dread(dataset_id,H5T_NATIVE_DOUBLE,memspace_id,dataspace_id,H5P_DEFAULT, hv);
00893
00894 #ifdef __DEBUG__VERBOSE__
00895 cout << "[DEBUGOUT] : " << __FILE__ << " : " << __LINE__ << endl;
00896 #endif
00899 hfield=new double*[nsample_];
00900
00901 for(i=0;i<nsample_;i++)
00902 {
00903 hfield[i]=new double[ncol];
00904
00905 for(j=0;j<ncol;j++)
00906 {
00907 hfield[i][j]=hv[i*ncol+j];
00908 }
00909 }
00910
00912 delete []hv;
00913
00915 H5Dclose(dataset_id);
00916 H5Sclose(dataspace_id);
00917 H5Sclose(memspace_id);
00918 H5Gclose(group_id);
00919 H5Fclose(file_id);
00920
00921
00922 return(OKCODE);
00923 }
00924
00925
00927 int h5_read_eigenvalues(std::string filename, std::string groupname,double* &lambdafreqquality)
00928 {
00930 int rank;
00931 int nmode=0;
00932 hsize_t dimsf[2];
00933 hid_t plist_id, file_id, group_id, dataset_id, dataspace_id, filespace_id, memspace_id;
00934 herr_t status;
00935
00936 plist_id = H5Pcreate(H5P_FILE_ACCESS);
00937 file_id = H5Fopen(filename.c_str(), H5F_ACC_RDONLY, plist_id);
00938 assert(file_id >= 0);
00939 H5Pclose(plist_id);
00940
00941 #ifdef HAVE_RLOG
00942 rInfo("...opening group in hdf5 file that stores the eigenvalues, eigenfrequencies and quality figures");
00943 #endif // HAVE_RLOG
00944
00945 group_id = H5Gopen(file_id, groupname.c_str());
00946
00947 #ifdef HAVE_RLOG
00948 rInfo("opened group in hdf5 file that stores the eigenvalues, eigenfrequencies and quality figures");
00949 #endif // HAVE_RLOG
00950
00952 dataset_id = H5Dopen(group_id,H5_NAME_FEMAX_EIGENVALUE.c_str());
00953 dataspace_id = H5Dget_space(dataset_id);
00954 rank = H5Sget_simple_extent_ndims(dataspace_id);
00955 status = H5Sget_simple_extent_dims(dataspace_id, dimsf,NULL);
00956 filespace_id = H5Screate_simple(2, dimsf, NULL);
00957 nmode=dimsf[0];
00958
00959 #ifdef HAVE_RLOG
00960 rInfo("eigenmodal solution contains %d modes", nmode);
00961 #endif // HAVE_RLOG
00962
00963
00965 hsize_t offset[2];
00966 hsize_t count[2];
00967
00968 offset[0] = 0;
00969 offset[1] = 0;
00970 count[0] = nmode;
00971 count[1] = dimsf[1];
00972
00974 status = H5Sselect_hyperslab(dataspace_id, H5S_SELECT_SET, offset, NULL, count, NULL);
00975
00977 hsize_t dimsm[2];
00978
00979 dimsm[0]=dimsf[0];
00980 dimsm[1]=dimsf[1];
00981 memspace_id=H5Screate_simple(rank,dimsm,NULL);
00982
00984 hsize_t offsetm[2];
00985 hsize_t countm[2];
00986
00987 offsetm[0]=0;
00988 offsetm[1]=0;
00989
00990 countm[0]=dimsf[0];
00991 countm[1]=dimsf[1];
00992
00993 status=H5Sselect_hyperslab(memspace_id,H5S_SELECT_SET, offsetm, NULL, countm, NULL);
00994
00996 lambdafreqquality=new double[nmode*3];
00997
00999 status=H5Dread(dataset_id,H5T_NATIVE_DOUBLE,memspace_id,dataspace_id,H5P_DEFAULT, lambdafreqquality);
01000
01002 H5Dclose(dataset_id);
01003 H5Sclose(dataspace_id);
01004 H5Sclose(memspace_id);
01005 H5Gclose(group_id);
01006 H5Fclose(file_id);
01007
01010
01011 return(nmode);
01012 }
01013
01015 int vtk_write_mesh(std::string vtklffname, int npoint, double** coord3d, int ntet, id_t** tetvertid)
01016 {
01017 int i,j;
01020 vtklffname.append(VTK_LFF_EXTENSION.c_str());
01021 std::ofstream of;
01022
01023 of.open(vtklffname.c_str());
01024 assert(of.is_open());
01025
01027 of.precision(DEFAULT_FLOAT_PRECISION);
01028
01030 of << "# vtk DataFile Version 2.0" << endl;
01031 of << "generated by femax_post3d" << endl;
01032 of << "ASCII" << endl << endl;
01033 of << "DATASET UNSTRUCTURED_GRID" << endl;
01034 of << "POINTS " << npoint << " float" << endl;
01035 of << scientific;
01036
01037 for (i = 0; i < npoint ; i++)
01038 {
01039 for(j=0;j<NREALCOORD3D;j++)
01040 {
01041 of << coord3d[i][j] << " ";
01042 }
01043 of << endl;
01044 }
01045
01046 of << endl;
01047 of << "CELLS " << ntet << SINGLE_SPACE.c_str() << 5*ntet << endl;
01048
01049
01050
01051 for (i = 0; i < ntet; i ++)
01052 {
01053 of << NCORNERTET;
01055 for (j = 0; j < NCORNERTET; j ++)
01056 {
01057 of << " " << tetvertid[i][j] << TAB_STRING.c_str();
01058 }
01059
01060 of << endl;
01061 }
01062
01063 of << endl;
01064
01065
01066 of << "CELL_TYPES " << ntet << endl;
01067 for ( i = 0; i < ntet; i++)
01068 {
01069 of << VTK_LFF_TETRAHEDRON_CODE << endl;
01070 }
01071
01072 of << endl;
01073
01074
01075 of << "CELL_DATA " << ntet << endl;
01076 of << "SCALARS process" << " float" << endl;
01077 of << "LOOKUP_TABLE default" << endl;
01078
01079 int processidx=NCORNERTET;
01081 for (i =0; i < ntet; i++)
01082 {
01083 of << tetvertid[i][processidx] << endl;
01084 }
01085 of << endl;
01086
01088 of.close();
01089
01090
01091 return(OKCODE);
01092 }
01093
01094
01096 int vtk_write_field(std::string vtklffname,
01097 int npoint,
01098 double** coord3d,
01099 int ntet,
01100 id_t** tetvertid,
01101 double*** efield,
01102 double*** hfield,
01103 int nmode)
01104 {
01105 int i,j,m;
01108 vtklffname.append(VTK_LFF_EXTENSION.c_str());
01109 std::ofstream of;
01110
01111 of.open(vtklffname.c_str());
01112 assert(of.is_open());
01113
01115 of.precision(DEFAULT_FLOAT_PRECISION);
01116
01119 #ifdef HAVE_RLOG
01120 rInfo("...writing VTK file header");
01121 #endif // HAVE_RLOG
01122
01123 of << "# vtk DataFile Version 2.0" << endl;
01124 of << "produce by femax_post3d - field and mesh solution" << endl;
01125 of << "ASCII" << endl << endl;
01126 of << "DATASET UNSTRUCTURED_GRID" << endl;
01127 of << "POINTS " << npoint << " float" << endl;
01128 of << scientific;
01129
01130 #ifdef HAVE_RLOG
01131 rInfo("wrote VTK file header");
01132 rInfo("...writing point information to VTK file");
01133 #endif // HAVE_RLOG
01134
01135 for (i = 0; i < npoint ; i++)
01136 {
01137 for(j=0;j<NREALCOORD3D;j++)
01138 {
01139 of << coord3d[i][j] << " ";
01140 }
01141 of << endl;
01142 }
01143
01144 of << endl;
01145
01146 #ifdef HAVE_RLOG
01147 rInfo("wrote VTK point information");
01148 rInfo("...writing VTK cell information");
01149 #endif // HAVE_RLOG
01150
01151 of << "CELLS " << ntet << SINGLE_SPACE.c_str() << 5*ntet << endl;
01152
01153
01154
01155 for (i = 0; i < ntet; i ++)
01156 {
01157 of << NCORNERTET;
01159 for (j = 0; j < NCORNERTET; j ++)
01160 {
01161 of << " " << tetvertid[i][j] << TAB_STRING.c_str();
01162 }
01163
01164 of << endl;
01165 }
01166
01167 of << endl;
01168
01169
01170 of << "CELL_TYPES " << ntet << endl;
01171 for ( i = 0; i < ntet; i++)
01172 {
01173 of << VTK_LFF_TETRAHEDRON_CODE << endl;
01174 }
01175 of << endl;
01176
01178 of << "CELL_DATA " << ntet << endl;
01179
01180 #ifdef HAVE_RLOG
01181 rInfo("wrote VTK CELL information");
01182 rInfo("...writing VTK VECTOR cell data for %d modes", nmode);
01183 #endif // HAVE_RLOG
01184
01186 for(m=0;m<nmode;m++)
01187 {
01188
01189 #ifdef HAVE_RLOG
01190 rInfo(" ");
01191 rInfo("...writing electric field VTK VECTOR cell data for mode %d", m);
01192 #endif // HAVE_RLOG
01193
01194
01195 of << "VECTORS efieldM" << m << " float" << endl;
01196
01197 for (i =0; i < ntet; i++)
01198 {
01199 for(j=0;j<NREALCOORD3D;j++)
01200 {
01201 of << efield[m][i][j] << TAB_STRING.c_str();
01202 }
01203 of << endl;
01204 }
01205 of << endl;
01206
01207 #ifdef HAVE_RLOG
01208 rInfo("wrote electric field VTK VECTOR cell data for mode %d", m);
01209 #endif // HAVE_RLOG
01210
01211 }
01212
01213
01214 of << endl;
01215
01216 #ifdef HAVE_RLOG
01217 rInfo(" ");
01218 #endif // HAVE_RLOG
01219
01220
01222 for(m=0;m<nmode;m++)
01223 {
01224 #ifdef HAVE_RLOG
01225 rInfo(" ");
01226 rInfo("...writing magnetic field VTK VECTOR cell data for mode %d", m);
01227 #endif // HAVE_RLOG
01228
01229 of << "VECTORS hfieldM" << m << " float" << endl;
01230
01231 for (i =0; i < ntet; i++)
01232 {
01233 for(j=0;j<NREALCOORD3D;j++)
01234 {
01235 of << hfield[m][i][j] << TAB_STRING.c_str();
01236 }
01237 of << endl;
01238 }
01239 of << endl;
01240
01241 #ifdef HAVE_RLOG
01242 rInfo("wrote magnetic field VTK VECTOR cell data for mode %d", m);
01243 #endif // HAVE_RLOG
01244
01245 }
01246
01247 #ifdef HAVE_RLOG
01248 rInfo(" ");
01249 #endif // HAVE_RLOG
01250
01251
01253 of.close();
01254
01255
01256 return(OKCODE);
01257 }
01258
01259
01261 int scf_write_efield(std::string scfname, int npoint, double** sloc, double*** efield, int nmode)
01262 {
01263 int i,j;
01266 scfname.append(".electric");
01267 scfname.append(SCF_EXTENSION.c_str());
01268 std::ofstream of;
01269
01270 of.open(scfname.c_str());
01271 assert(of.is_open());
01272
01274 of.precision(DEFAULT_FLOAT_PRECISION);
01275
01278 #ifdef HAVE_RLOG
01279 rInfo("...writing header of ANSYS readable ASCII file storing the electric field");
01280 #endif // HAVE_RLOG
01281
01282
01283 of << "# simple column format for use by Markus Bopp - Version 2006 " << endl;
01284 of << "# generated by femax_post3d" << endl;
01285 of << "# format: <x> <y> <z> <ex0> <ey0> <ez0> <ex0> <ey0> <ez0> ... <ex N-1> <ey N-1> <ez N-1> where N = number of computed modes" << endl;
01286 of << "# nsamples " << npoint << endl;
01287 of << scientific;
01288
01289 #ifdef HAVE_RLOG
01290 rInfo("wrote header of ANSYS readable ASCII file for the electric field");
01291 #endif // HAVE_RLOG
01292
01293
01294
01295
01296 #ifdef HAVE_RLOG
01297 rInfo("...writing electric field for all modes into ANSYS readable ASCII file");
01298 #endif // HAVE_RLOG
01299
01300
01301 for (i = 0; i < npoint ; i++)
01302 {
01303
01304 for(j=0;j<NREALCOORD3D;j++)
01305 {
01306 of << sloc[i][j] << " ";
01307 }
01308
01309 of << TAB_STRING.c_str();
01310
01311
01312 for(unsigned int m=0;m<nmode;m++)
01313 {
01314 for(j=0;j<NREALCOORD3D;j++)
01315 {
01316 of << efield[m][i][j] << TAB_STRING.c_str();
01317 }
01318 }
01319
01320
01321 of << endl;
01322 }
01323
01324 #ifdef HAVE_RLOG
01325 rInfo("wrote electric field for all modes into ANSYS readable ASCII file");
01326 #endif // HAVE_RLOG
01327
01328
01329
01331 of.close();
01332
01333
01334 return(OKCODE);
01335 }
01336
01338 int scf_write_hfield(std::string scfname, int npoint, double** sloc, double*** hfield, int nmode)
01339 {
01340 int i,j,m;
01343 scfname.append(".magnetic");
01344 scfname.append(SCF_EXTENSION.c_str());
01345 std::ofstream of;
01346
01347 of.open(scfname.c_str());
01348 assert(of.is_open());
01349
01351 of.precision(DEFAULT_FLOAT_PRECISION);
01352
01355 #ifdef HAVE_RLOG
01356 rInfo("...writing header of ANSYS readable ASCII file storing the magnetic field");
01357 #endif // HAVE_RLOG
01358
01359 of << "# simple column format for use by Markus Bopp - Version 2006 " << endl;
01360 of << "# generated by femax_post3d" << endl;
01361 of << "# format: <x> <y> <z> <ex0> <ey0> <ez0> <ex0> <ey0> <ez0> ... <ex N-1> <ey N-1> <ez N-1> where N = number of computed modes" << endl;
01362 of << "# nsamples " << npoint << endl;
01363 of << scientific;
01364
01365 #ifdef HAVE_RLOG
01366 rInfo("wrote header of ANSYS readable ASCII file for the magnetic field");
01367 #endif // HAVE_RLOG
01368
01369
01370
01371
01372 #ifdef HAVE_RLOG
01373 rInfo("...writing magnetic field for all modes into ANSYS readable ASCII file");
01374 #endif // HAVE_RLOG
01375
01376 for (i = 0; i < npoint ; i++)
01377 {
01378
01379 for(j=0;j<NREALCOORD3D;j++)
01380 {
01381 of << sloc[i][j] << " ";
01382 }
01383
01384 of << TAB_STRING.c_str();
01385
01386
01387 for(m=0;m<nmode;m++)
01388 {
01389 for(j=0;j<NREALCOORD3D;j++)
01390 {
01391 of << hfield[m][i][j] << TAB_STRING.c_str();
01392 }
01393 }
01394
01395
01396 of << endl;
01397 }
01398
01399 #ifdef HAVE_RLOG
01400 rInfo("wrote magnetic field for all modes into ANSYS readable ASCII file");
01401 #endif // HAVE_RLOG
01402
01404 of.close();
01405
01406
01407 return(OKCODE);
01408 }
01409
01410
01411
01412