advanced/dsbandrepair/dsbandrepair.cc

0001 //
0002 // ********************************************************************
0003 // * License and Disclaimer                                           *
0004 // *                                                                  *
0005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
0006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
0007 // * conditions of the Geant4 Software License,  included in the file *
0008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
0009 // * include a list of copyright holders.                             *
0010 // *                                                                  *
0011 // * Neither the authors of this software system, nor their employing *
0012 // * institutes,nor the agencies providing financial support for this *
0013 // * work  make  any representation or  warranty, express or implied, *
0014 // * regarding  this  software system or assume any liability for its *
0015 // * use.  Please see the license in the file  LICENSE  and URL above *
0016 // * for the full disclaimer and the limitation of liability.         *
0017 // *                                                                  *
0018 // * This  code  implementation is the result of  the  scientific and *
0019 // * technical work of the GEANT4 collaboration.                      *
0020 // * By using,  copying,  modifying or  distributing the software (or *
0021 // * any work based  on the software)  you  agree  to acknowledge its *
0022 // * use  in  resulting  scientific  publications,  and indicate your *
0023 // * acceptance of all terms of the Geant4 Software license.          *
0024 // ********************************************************************
0025 //
0026 //  author: Le Tuan Anh, 20/10/2023
0027 /// \file main.cc
0028 /// \brief Main program of the dsbandrepair
0029
0030 #include "G4UImanager.hh"
0031 #include "G4UIterminal.hh"
0032 #include "G4UItcsh.hh"
0033 #include "G4UIExecutive.hh"
0034
0035 #include "G4RunManagerFactory.hh"
0036
0037 #ifdef G4VIS_USE
0038 #include "G4VisExecutive.hh"
0039 #endif
0040
0041 #include "G4Timer.hh"
0042 #include "G4ExceptionSeverity.hh"
0043 #include "G4DNAChemistryManager.hh"
0044 #include "G4VisExecutive.hh"
0045 #include "G4Filesystem.hh"
0046
0047 #include "PhysActionInitialization.hh"
0048 #include "DetectorConstruction.hh"
0049 #include "PhysicsList.hh"
0050 #include "InformationKeeper.hh"
0051
0052 #include "ChemActionInitialization.hh"
0053 #include "ChemPhysicsList.hh"
0054 #include "ChemNtupleManager.hh"
0055
0056 #ifdef USE_MPI
0057 #include "G4MPImanager.hh"
0058 #include "G4MPIsession.hh"
0059 #include "G4MPIextraWorker.hh"
0060 #endif
0061 #include <ctime>
0062
0063 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
0064
0065 void CheckingSomeFilesAndFolders();
0066 G4String ExtractChemListNameFromMacroFile(G4String);
0067
0068 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
0069
0070 int main(int argc,char** argv)
0071 {
0072 #ifdef USE_MPI
0073         G4MPImanager* g4MPI = new G4MPImanager(argc, argv, 0);
0074         g4MPI->SetVerbose(1);
0075         G4MPIsession* session = g4MPI-> GetMPIsession();
0076         G4String prompt = "[40;01;33m";
0077         prompt += "G4MPI";
0078         prompt += "[40;31m(%s)[40;36m[%/][00;30m:";
0079         session-> SetPrompt(prompt);
0080 #else
0081         G4UIExecutive* ui = nullptr;
0082         if ( argc == 1 ) { ui = new G4UIExecutive(argc, argv); }
0083 #endif // USE_MPI
0084     if (argc < 2) {
0085         G4cerr<<"====>> Wrong input. To run Physgeo, type : ./dsbandrepair macrofile\n"
0086             <<"To run Chem_geo, type : ./dsbandrepair macrofile chem"<<G4endl;
0087 #ifdef USE_MPI
0088         delete g4MPI;
0089 #endif // USE_MPI
0090         return EXIT_FAILURE;
0091     }
0092     G4String macrofileName = argv[1];
0093     if (argc > 2) {
0094         const G4String rmode = argv[2];
0095         if (rmode == "phys") gRunMode = RunningMode::Phys;
0096         if (rmode == "chem") gRunMode = RunningMode::Chem;
0097     }
0098     // Choose the Random engine
0099     time_t timeStart;
0100     time(&timeStart);
0101     unsigned long seed = timeStart;
0102 #ifdef USE_MPI
0103     // Le Tuan Anh: add rankID to get different seeds for multi parallel processes
0104     seed += 1987*g4MPI->GetRank();
0105 #endif // USE_MPI
0106     G4cout<<"Initial Seed for random engine: "<<seed<<G4endl;
0107     CLHEP::HepRandom::setTheEngine(new CLHEP::MTwistEngine);
0108     CLHEP::HepRandom::setTheSeed(seed);
0109     if (gRunMode == RunningMode::Phys) {
0110 #ifdef USE_MPI
0111         if (g4MPI->GetRank() == 0 ){
0112             CheckingSomeFilesAndFolders();
0113         }
0114 #else
0115         CheckingSomeFilesAndFolders();
0116 #endif // USE_MPI
0117         G4Random::setTheSeed(seed);
0118         auto* runManager = G4RunManagerFactory::CreateRunManager(G4RunManagerType::Default);
0119 #ifdef G4MULTITHREADED
0120         G4int threadNumber= 1;
0121         runManager-> SetNumberOfThreads(threadNumber);
0122 #endif  // G4MULTITHREADED
0123
0124         DetectorConstruction* detector = new DetectorConstruction(1.,0,false);
0125         runManager->SetUserInitialization(detector);
0126
0127         PhysicsList* physList = new PhysicsList;
0128         runManager->SetUserInitialization(physList);
0129
0130         PhysActionInitialization* actionIni = new PhysActionInitialization();
0131         runManager->SetUserInitialization(actionIni);
0132 #ifdef USE_MPI
0133         // extra worker (for collecting ntuple data)
0134         if ( g4MPI->IsExtraWorker() ) {
0135             G4cout << "Set extra worker" << G4endl;
0136             G4UserRunAction* runAction = const_cast<G4UserRunAction*>(runManager->GetUserRunAction());
0137             g4MPI->SetExtraWorker(new G4MPIextraWorker(runAction));
0138         }
0139         session-> SessionStart();
0140
0141
0142         if (g4MPI->GetRank() == 0 ){
0143             InformationKeeper::Instance()->WritePhysGeo();
0144         }
0145         delete g4MPI;
0146 #else
0147         // Get the pointer to the User Interface manager
0148         G4UImanager* UImanager = G4UImanager::GetUIpointer();
0149         // Process macro or start UI session
0150         if ( ! ui ) {
0151             // batch mode
0152             G4String command = "/control/execute ";
0153             UImanager->ApplyCommand(command+macrofileName);
0154         }
0155         InformationKeeper::Instance()->WritePhysGeo();
0156 #endif // USE_MPI
0157         delete runManager;
0158     } else if (gRunMode == RunningMode::Chem) {
0159         std::string inputFileorFolder = "chem_input";
0160         if (argc == 4) inputFileorFolder = argv[3];
0161         G4fs::path p{inputFileorFolder};
0162         G4String outputFileName = "test";
0163         std::vector<G4String> totalNumberofFilesVector, numberOfFilesTobeProcessedVector;
0164         if (G4fs::is_directory(p)) {
0165             for (const auto& entry : G4fs::directory_iterator(p)) {
0166                 if (entry.path().extension() == ".dat") {
0167                     totalNumberofFilesVector.push_back(entry.path().string());
0168                 }
0169             }
0170             std::sort(totalNumberofFilesVector.begin(),totalNumberofFilesVector.end());
0171 #ifdef USE_MPI
0172             G4int numberofRanks = g4MPI->GetActiveSize();
0173             size_t filesTobeProcessedSlave = (size_t)(
0174                 std::floor(G4double(totalNumberofFilesVector.size())/G4double(numberofRanks)));
0175             // note: should not use "std::ceil"
0176             size_t filesTobeProcessedMaster =
0177                 totalNumberofFilesVector.size() - (numberofRanks-1)*filesTobeProcessedSlave;
0178             if (g4MPI->IsMaster()) {
0179                 for (size_t ii=0; ii< filesTobeProcessedMaster; ii++) {
0180                     numberOfFilesTobeProcessedVector.push_back(totalNumberofFilesVector.at(ii));
0181                 }
0182             } else {
0183                 for (size_t ii=0; ii< filesTobeProcessedSlave; ii++) {
0184                     auto rankID = g4MPI->GetRank();
0185                     size_t kk = filesTobeProcessedMaster + (rankID-1)*filesTobeProcessedSlave + ii;
0186                     numberOfFilesTobeProcessedVector.push_back(totalNumberofFilesVector.at(kk));
0187                 }
0188             }
0189             G4cout<<"-----> "<<numberOfFilesTobeProcessedVector.size()
0190                     <<" files will be processed on rank #"<<g4MPI->GetRank()<<G4endl;
0191 #else
0192             numberOfFilesTobeProcessedVector = totalNumberofFilesVector;
0193 #endif
0194             if (totalNumberofFilesVector.size() == 0) {
0195                 G4cout<<"===>> There is no files found in "<<inputFileorFolder
0196                 <<". You have to run Phys_geo first!!!"<<G4endl;
0197                 //return EXIT_SUCCESS;
0198             } else {
0199                 G4cout<<"===>> Total files found in "<<inputFileorFolder
0200                 <<" : "<<totalNumberofFilesVector.size()<<G4endl;
0201             }
0202         } else if (G4fs::is_regular_file(p)) {
0203             numberOfFilesTobeProcessedVector.push_back(inputFileorFolder);
0204             if (p.has_stem()) {
0205                 outputFileName = p.stem().string();
0206             } else outputFileName = inputFileorFolder;
0207         }
0208         else G4cout<<"===>>dsbandrepair: "<<p.string()<<" is Not Directory or file !!!"<<G4endl;
0209
0210         //------------------------------------------
0211         // Initialization classes
0212         //------------------------------------------
0213
0214         auto* runManager = G4RunManagerFactory::CreateRunManager(G4RunManagerType::Serial);
0215
0216         G4DNAChemistryManager::Instance()->SetChemistryActivation(true);
0217
0218         ChemNtupleManager ntupleManager;
0219
0220         DetectorConstruction* detector = new DetectorConstruction();
0221         G4String firstFileForInit="";
0222         if (numberOfFilesTobeProcessedVector.size()>0) {
0223             firstFileForInit=numberOfFilesTobeProcessedVector.at(0);
0224             detector->ParseGeoFileForChemMode(firstFileForInit); // read to build voxel
0225         }
0226         runManager->SetUserInitialization(detector);
0227         G4String chemListName = ExtractChemListNameFromMacroFile(macrofileName);
0228         ChemPhysicsList* physList;
0229         if (chemListName != "") {
0230             physList = new ChemPhysicsList(chemListName);
0231         } else physList = new ChemPhysicsList();
0232         runManager->SetUserInitialization(physList);
0233
0234         ChemActionInitialization* actionIni = new ChemActionInitialization(&ntupleManager,physList);
0235         runManager->SetUserInitialization(actionIni);
0236
0237         //get the pointer to the User Interface manager
0238         G4UImanager* UI = G4UImanager::GetUIpointer();
0239         G4String command = "/control/execute ";
0240         UI->ApplyCommand(command+macrofileName);
0241         runManager->Initialize();
0242         if (numberOfFilesTobeProcessedVector.size()>0) {
0243             size_t nprocessedfiles{0}, ncounts{1};
0244             if (numberOfFilesTobeProcessedVector.size()>=100) ncounts=10;
0245             if (numberOfFilesTobeProcessedVector.size()>=10000) ncounts=100;
0246             if (numberOfFilesTobeProcessedVector.size()>=1000000) ncounts=500;
0247             for (auto const &fileInput : numberOfFilesTobeProcessedVector) {
0248                 G4fs::path aP{std::string(fileInput)};
0249                 if (aP.has_stem()) {
0250                     outputFileName = aP.stem().string();
0251                 } else outputFileName = fileInput;
0252                 ntupleManager.SetFileName(outputFileName);
0253                 if (fileInput != firstFileForInit) detector->ParseGeoFileForChemMode(fileInput);
0254                 detector->InsertMoleculeInWorld();
0255                 UI->ApplyCommand("/run/beamOn 1");
0256                 nprocessedfiles++;
0257                 if (nprocessedfiles == 1 ||
0258                     nprocessedfiles == numberOfFilesTobeProcessedVector.size() ||
0259                     0 == (nprocessedfiles % ncounts)) {
0260                     G4cout<<"=====> Processed file: "<<nprocessedfiles<<"-th/("
0261                           <<numberOfFilesTobeProcessedVector.size()<<" files)"
0262 #ifdef USE_MPI
0263                           <<" in rank #"<<g4MPI->GetRank()
0264 #endif
0265                           <<"!!!"<<G4endl;
0266                 }
0267             }
0268         } else {
0269             UI->ApplyCommand("/run/beamOn 1");
0270         }
0271 #ifdef USE_MPI
0272         delete g4MPI;
0273 #endif // USE_MPI
0274         delete runManager;
0275     } else {
0276         G4cout<<"Undefined Running Mode; dsbansrepair will quit now. See you!\n";
0277     }
0278
0279     return EXIT_SUCCESS;
0280 }
0281
0282 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
0283
0284 void CheckingSomeFilesAndFolders()
0285 {
0286     const G4fs::path curPath{G4fs::current_path()};
0287     // create output folder for containing Phys output
0288     G4fs::file_status fst = G4fs::file_status{};
0289     const G4fs::path outFolderPhysP{InformationKeeper::Instance()->GetPhysOutFolderName().c_str()};
0290     auto isExist = G4fs::status_known(fst) ? G4fs::exists(fst) : G4fs::exists(outFolderPhysP);
0291     if (! isExist) {
0292         G4fs::create_directory(outFolderPhysP);
0293     } else {
0294         G4fs::remove_all(outFolderPhysP);//delete old folder
0295         G4fs::create_directory(outFolderPhysP);
0296     }
0297     // create output folder for containing chem_input
0298     const G4fs::path outFolderP{InformationKeeper::Instance()->GetChemInputFolderName().c_str()};
0299     isExist = G4fs::status_known(fst) ? G4fs::exists(fst) : G4fs::exists(outFolderP);
0300     if (! isExist) {
0301         G4fs::create_directory(outFolderP);
0302     } else {
0303         G4fs::remove_all(outFolderP);//delete old folder
0304         G4fs::create_directory(outFolderP);
0305     }
0306 }
0307
0308 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0309
0310 G4String ExtractChemListNameFromMacroFile(G4String fileName)
0311 {
0312     G4String out="";
0313     std::ifstream file;
0314     file.open(fileName.c_str());
0315     if (!file.is_open()) {
0316         G4String msg = "Error in openning file " + fileName;
0317         G4Exception("G4String ExtractChemListNameFromMacroFile()", "", FatalException, msg);
0318     } else {
0319         G4String line;
0320         while(std::getline(file, line))
0321         {
0322             std::istringstream iss(line);
0323             G4String flag;
0324             iss >> flag;
0325             G4String tvalue;
0326             iss >> tvalue;
0327             if (flag == "/dsbandrepair/chem/chemList") out = tvalue;
0328         }
0329     }
0330     file.close();
0331     return out;
0332 }
0333
0334 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......