!*********************************************************************** ! ! Copyright (c) 2015, Lawrence Livermore National Security, LLC. ! Produced at the Lawrence Livermore National ! Laboratory. ! Written by Nicolas Schunck, schunck1@llnl.gov ! ! LLNL-CODE-470611 All rights reserved. ! ! Copyright 2011, N. Schunck, J. Dobaczewski, J. McDonnell, ! W. Satula, J.A. Sheikh, A. Staszczak, ! M. Stoitsov, P. Toivanen ! Copyright 2009, J. Dobaczewski, W. Satula, B.G. Carlsson, J. Engel, ! P. Olbratowski, P. Powalowski, M. Sadziak, ! J. Sarich, N. Schunck, A. Staszczak, M. Stoitsov, ! M. Zalewski, H. Zdunczuk ! Copyright 2004, 2005, J. Dobaczewski, P. Olbratowski ! Copyright 1997, 2000, J. Dobaczewski, J. Dudek ! ! This file is part of HFODD. ! ! HFODD is free software: you can redistribute it and/or modify it ! under the terms of the GNU General Public License as published by ! the Free Software Foundation, either version 3 of the License, or ! (at your option) any later version. ! ! HFODD is distributed in the hope that it will be useful, but ! WITHOUT ANY WARRANTY; without even the implied warranty of ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the ! GNU General Public License for more details. ! ! You should have received a copy of the GNU General Public License ! along with HFODD. If not, see . ! ! OUR NOTICE AND TERMS AND CONDITIONS OF THE GNU GENERAL PUBLIC ! LICENSE ! ! Our Preamble Notice ! ! A. This notice is required to be provided under our contract ! with the U.S. Department of Energy (DOE). This work was ! produced at the Lawrence Livermore National Laboratory under ! Contract No. DE-AC52-07NA27344 with the DOE. ! B. Neither the United States Government nor Lawrence Livermore ! National Security, LLC nor any of their employees, makes any ! warranty, express or implied, or assumes any liability or ! responsibility for the accuracy, completeness, or usefulness ! of any information, apparatus, product, or process disclosed, ! or represents that its use would not infringe privately-owned ! rights. ! C. Also, reference herein to any specific commercial products, ! process, or services by trade name, trademark, manufacturer ! or otherwise does not necessarily constitute or imply its ! endorsement, recommendation, or favoring by the United States ! Government or Lawrence Livermore National Security, LLC. The ! views and opinions of authors expressed herein do not ! necessarily state or reflect those of the United States ! Government or Lawrence Livermore National Security, LLC, and ! shall not be used for advertising or product endorsement ! purposes. ! ! The precise terms and conditions for copying, distribution and ! modification are contained in the file COPYING. ! !*********************************************************************** !----------------------------------------------------------------------! ! ! ! MPI I/O Environment for HFODD ! ! ! ! Authors: N. Schunck, ORNL ! ! ! ! This module is/will be in charge of: ! ! - reading sequential and parallel input data from the hfodd.d ! ! and hfodd_mpiio.d files ! ! - broadcasting all input data to all processors ! ! - converting data to the proper variables used in HFODD ! ! ! ! First included in official release v249e ! ! ! !----------------------------------------------------------------------! Module linkedLists Implicit None ! Define derived types Node_int and Node_dbl, which are used ! to construct linked list of integers and real numbers Type :: Node_int Integer :: value Integer :: numero Character(Len=10) :: name TYPE(Node_int), pointer :: next End type Node_int Type :: Node_dbl Real :: value Integer :: numero Character(Len=10) :: name TYPE(Node_dbl), pointer :: next End type Node_dbl Type :: Node_str Character(Len=10) :: value Integer :: numero Character(Len=10) :: name TYPE(Node_str), pointer :: next End type Node_str Contains ! A bunch of useful subroutines to manipulate linked lists ! Routines dealing with Node_int type Recursive Subroutine printEntireList_int(ptr) Type(Node_int), pointer :: ptr If(associated(ptr%next)) Call printEntireList_int(ptr%next) print *, ptr%value, ptr%numero, ptr%name End Subroutine printEntireList_int Recursive Subroutine freeEntireList_int(ptr) Type(Node_int), pointer :: ptr If(associated(ptr%next)) Call freeEntireList_int(ptr%next) Deallocate(ptr) End Subroutine freeEntireList_int ! Routines dealing with Node_dbl type Recursive Subroutine printEntireList_dbl(ptr) Type(Node_dbl), pointer :: ptr If(associated(ptr%next)) Call printEntireList_dbl(ptr%next) Print *, ptr%value, ptr%numero, ptr%name End Subroutine printEntireList_dbl Recursive Subroutine freeEntireList_dbl(ptr) Type(Node_dbl), pointer :: ptr If(associated(ptr%next)) Call freeEntireList_dbl(ptr%next) Deallocate(ptr) End Subroutine freeEntireList_dbl ! Routines dealing with Node_str type Recursive Subroutine printEntireList_str(ptr) Type(Node_str), pointer :: ptr If(associated(ptr%next)) Call printEntireList_str(ptr%next) Print *, ptr%value, ptr%numero, ptr%name End Subroutine printEntireList_str Recursive Subroutine freeEntireList_str(ptr) Type(Node_str), pointer :: ptr If(associated(ptr%next)) Call freeEntireList_str(ptr%next) Deallocate(ptr) End Subroutine freeEntireList_str End Module linkedLists Module hfodd_mpiio Use linkedLists Implicit None Include "mpif.h" ! Private Integer, PRIVATE :: stringLength = 10 Integer, PRIVATE :: counterInt, counterDbl, counterStr, counterTot Integer, PRIVATE :: counterVector(1:3) Integer, Allocatable, PRIVATE :: vectorInt_value(:) Real, Allocatable, PRIVATE :: vectorDbl_value(:) Character(Len=10), Allocatable, PRIVATE :: vectorStr_value(:) Character(Len=10), Allocatable, PRIVATE :: vectorAll_names(:) Type(Node_int), pointer :: listeInt Type(Node_dbl), pointer :: listeDbl Type(Node_str), pointer :: listeStr ! Public Integer, PUBLIC :: modeMPI_deformation,modeMPI_neck,modeMPI_basis Integer, PUBLIC :: numberConstraints,numberConstraints_restart,restartTemperature, & neck_numberPoints,neck_numberPoints_restart,counterForce, & counterTemperature,batch_mode,number_batch,number_constraint_batch Integer, PUBLIC :: neutron_stepNumber,neutron_stepValue,neutron_initValue, & proton_stepNumber,proton_stepValue,proton_initValue, & basis_frequency_stepNumber,basis_deform_stepNumber, & basis_size_stepNumber,basis_size_stepValue,basis_size_initValue, & basis_states_stepNumber,basis_states_stepValue,basis_states_initValue, & turnConstraintOn,optimizePath Integer, Allocatable, PUBLIC :: qLambda(:),qMiu(:),qStepNumber(:),qLambda_restart(:), & qMiu_restart(:),qStepNumber_restart(:),lambda_batch(:),& miu_batch(:) Real, PUBLIC :: x_neutron_initValue,x_neutron_stepValue, & x_proton_initValue,x_proton_stepValue Real, PUBLIC :: basis_frequency_stepValue,basis_frequency_initValue, & basis_deform_stepValue,basis_deform_initValue, & neck_initValue,neck_finValue,neck_initValue_restart, & neck_finValue_restart Real, Allocatable, PUBLIC :: qinitValue(:),qFinValue(:),qinitValue_restart(:), & qFinValue_restart(:),temperatureVector(:) Character(Len=4), Allocatable, PUBLIC :: forceVector(:) Contains !------------------------------------------------------------------! ! This subroutine reads the standard HFODD input file and puts ! ! all encountered data into two linked lists of type Node_int ! ! (for data of type integer) and Node_dbl (for data of type real). ! ! These derived types also contain a character string that is ! ! used to identify the keyword where this data come from. This ! ! will serve for the reconstruction of the original variables. ! !------------------------------------------------------------------! Subroutine mpi_getSequentialData() Integer :: NFIREA Integer, Parameter :: NDITER = 5000 Integer :: ii, jj, counterDef, counterWS, counterColl Integer :: mpi_rank, mpi_err, mpi_string_10 Integer :: ICOTYP, ICOUDI, ICOUEX, I_GOGA, NOITER, IBROYD, N_DUMM, & NO_RUN, NULAST, NUPING, NUCHAO, ISIMPY, ISIGNY, IPARTY, & IROTAT, ISIMTX, ISIMTY, ISIMTZ, IPAIRI, IPAHFB, IPABCS, & IMFHFB, LIPKIN, LIPKIP, IOPTGS, NXHERM, NYHERM, NZHERM, & NOSCIL, NLIMIT, I1LINE, NILXYZ, IREVIE, ICONTI, IPCONT, & IGCONT, IGPCON, ILCONT, IFCONT, IACONT, IFSHEL, NPOLYN, & LAMBDA, ITE_WS, IREAWS, ISTAND, KETA_J, KETA_W, KETACM, & KETA_M, IF_RPA, IF_EDF, IF_THO, IBASIS, IHFPOT, IFRAGM, & L_PUSH, M_PUSH, N_PUSH, DOPUSH, NPLAST, MIU, NMUCON, & NMUCOU, NMUPRI, NEXBET, IPRIBE, IPRIBL, IWRIFI, INNUMB, & IZNUMB, IGOGPA, ITRMAX, NTHETA, MIN_QP, IDEALL, IDELOC, & IDECON, IF2FIN, IF2FIP Real :: A_DUMM, B_DUMM, ENECUT, SLOWEV, SLOWOD, SLOWPA, & SLOWLI, EPSITE, EPSPNG, ECUTOF, PRHO_T, PRHOST, & PRHO_N, PRHOSN, POWERN, PRHO_P, PRHOSP, POWERP, & POWERT, OMEGAY, GSTRUN, GSTRUP, HOMFAC, ALPHAR, & STP_WS, CBETHO, TEMP_T, Q_PUSH, FCHOM0, R0PARM, & DELTAE, XLOCMX, V2_MIN, FE2FIN, FE2FIP Real :: PI, DEGRAD Character(Len=10) :: KEYWOR Type(Node_int), Pointer :: tempNodeInt, currentNodeInt, previousNodeInt Type(Node_dbl), Pointer :: tempNodeDbl, currentNodeDbl, previousNodeDbl Type(Node_str), Pointer :: tempNodeStr, currentNodeStr, previousNodeStr ! Initialize local variables ICOTYP=7 ICOUDI=1;ICOUEX=1 I_GOGA=0;IGOGPA=0 NOITER=10 IBROYD=1 N_DUMM=8 NO_RUN=0 NULAST=3;NUPING=3;NUCHAO=3 ISIMPY=1;ISIGNY=1;IPARTY=1 IROTAT=0 ISIMTX=-1;ISIMTY=-1;ISIMTZ=-1 IPAIRI=1;IPAHFB=1 IPABCS=0 IMFHFB=0;LIPKIN=0;LIPKIP=0;IF2FIN=0;IF2FIP=0 IOPTGS=1 NXHERM=30;NYHERM=30;NZHERM=30 NOSCIL=12 NLIMIT=455 INNUMB=86;IZNUMB=66 LAMBDA=0; MIU=0 I1LINE=1 NILXYZ=3 IREVIE=0;ICONTI=0;IPCONT=0;IGCONT=0;IGPCON=0;ILCONT=0;IFCONT=0;IACONT=0 IFSHEL=0;NPOLYN=6 ITE_WS=10;IREAWS=0 ISTAND=1;KETA_J=0;KETA_W=1;KETACM=1;KETA_M=0 IF_RPA=1;IF_EDF=0;IF_THO=0;CBETHO=0.0D0 IBASIS=0 IFRAGM=0 ITRMAX=31;NTHETA=181 IDEALL=0;IDELOC=0;IDECON=0 IHFPOT=0 L_PUSH=1;M_PUSH=0;N_PUSH=0;DOPUSH=0;NPLAST=1 NMUCON=4;NMUCOU=4;NMUPRI=4 NEXBET=4;IPRIBE=1;IPRIBL=0 IWRIFI=0 A_DUMM=0.0D0; B_DUMM=0.0D0 ENECUT=800.0D0 R0PARM=1.2D0 ALPHAR=0.0D0 SLOWEV=0.5D0; SLOWOD=0.5D0; SLOWPA=0.5D0; SLOWLI=0.5D0 EPSITE=1.0D-5;EPSPNG=1.0d-5 ECUTOF=60.0D0 PRHO_T=-000.0D0; PRHOST=0.32D0; POWERT=1.0D0 PRHO_N=-200.0D0; PRHOSN=0.32D0; POWERN=1.0D0 PRHO_P=-200.0D0; PRHOSP=0.32D0; POWERP=1.0D0 OMEGAY=0.0D0 GSTRUN=1.2; GSTRUP=1.2;HOMFAC=3.5D0 STP_WS=0.01D0 Q_PUSH=0.00D0 FCHOM0=1.20D0 MIN_QP=0;DELTAE=2.50D0;XLOCMX=0.75D0;V2_MIN=0.005D0 FE2FIN=0.1D0;FE2FIP=0.1D0 ! Getting the rank of the current process Call mpi_comm_rank(MPI_COMM_WORLD, mpi_rank, mpi_err) ! Define MPI type for a character string of length stringLength and commit Call mpi_type_contiguous(stringLength, MPI_CHARACTER, mpi_string_10, mpi_err) Call mpi_type_commit(mpi_string_10, mpi_err) NO_RUN = 0 ! Only the root process reads the data If(mpi_rank.Eq.0) Then PI=4.0D0*Atan(1.0D0) DEGRAD=180.0D0/PI counterInt = 0 counterDbl = 0 counterStr = 0 ! Initializing the end of the linked list Write(KEYWOR,'(10x)') counterInt = counterInt + 1 Allocate(listeInt) listeInt%value = 0 listeInt%numero = counterInt listeInt%name = KEYWOR NULLIFY(listeInt%next) counterDbl = counterDbl + 1 Allocate(listeDbl) listeDbl%value = 0.0D0 listeDbl%numero = counterDbl listeDbl%name = KEYWOR NULLIFY(listeDbl%next) counterStr = counterStr + 1 Allocate(listeStr) listeStr%value = "INITIALIZE" listeStr%numero = counterStr listeStr%name = KEYWOR NULLIFY(listeStr%next) NFIREA=87 Open(NFIREA,FILE='hfodd.d',STATUS='UNKNOWN',FORM='FORMATTED') ! Main loop over the keyword value Do While (trim(KEYWOR) .NE. 'ALL_DONE' .AND. trim(KEYWOR) .NE. 'ALL-DONE') Read(NFIREA,'(A10)') KEYWOR !==============================! ! ITERATIONS CONTROL FLAGS ! !==============================! ! Number of iterations If(trim(KEYWOR).Eq.'ITERATIONS') Then Read(NFIREA,*) NOITER counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = NOITER tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt If(NOITER.GT.NDITER) STOP 'NOITER.GT.NDITER IN NAMELI' End If ! Switches on/off the Broyden method If(trim(KEYWOR).Eq.'BROYDEN') Then Read(NFIREA,*) IBROYD,N_DUMM,A_DUMM,B_DUMM counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IBROYD tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = N_DUMM tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = A_DUMM tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = B_DUMM tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl End If ! Slow-down parameters for the mean-field If(trim(KEYWOR).Eq.'SLOW_DOWN'.OR.trim(KEYWOR).Eq.'SLOW-DOWN') Then Read(NFIREA,*) SLOWEV,SLOWOD counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = SLOWEV tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = SLOWOD tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl End If ! Slow-down parameters for the pairing field If(trim(KEYWOR).Eq.'SLOW_PAIR'.OR.trim(KEYWOR).Eq.'SLOW-PAIR') Then Read(NFIREA,*) SLOWPA counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = SLOWPA tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl End If ! Slow-down parameters for the Lipkin-Nogami correction If(trim(KEYWOR).Eq.'SLOWLIPKIN') Then Read(NFIREA,*) SLOWLI counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = SLOWLI tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl End If ! Defines the precision of convergence If(trim(KEYWOR).Eq.'ITERAT_EPS'.OR.trim(KEYWOR).Eq.'ITERAT-EPS') Then Read(NFIREA,*) EPSITE counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = EPSITE tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl End If ! Defines divergence If(trim(KEYWOR).Eq.'MAXANTIOSC') Then Read(NFIREA,*) NULAST counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = NULAST tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Control ping-pong effects If(trim(KEYWOR).Eq.'PING_PONG'.OR.trim(KEYWOR).Eq.'PING-PONG') Then Read(NFIREA,*) EPSPNG,NUPING counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = NUPING tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = EPSPNG tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl End If ! Control chaotic divergence If(trim(KEYWOR).Eq.'CHAOTIC') Then Read(NFIREA,*) NUCHAO counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = NUCHAO tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If !=====================================! ! SPECIFIC FEATURES CONTROL FLAGS ! !=====================================! ! Definition of the collective space for calculation of collective inertia If(trim(KEYWOR).Eq.'COLL_SPACE'.OR.trim(KEYWOR).Eq.'COLL-SPACE') Then counterColl = 0 LAMBDA=-1 Do While (LAMBDA < 0) Read(NFIREA,*) LAMBDA,MIU counterColl = counterColl + 1 counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = ABS(LAMBDA) tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = MIU tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End Do counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = counterColl tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Type of shell correction If(trim(KEYWOR).Eq.'SHELLCORCT') Then Read(NFIREA,*) IFSHEL counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IFSHEL tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Shell correction parameters If(trim(KEYWOR).Eq.'SHELLPARAM') Then Read(NFIREA,*) GSTRUN,GSTRUP,HOMFAC,NPOLYN counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = GSTRUN tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = GSTRUP tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = HOMFAC tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = NPOLYN tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Switching to interface with HFBTHO If(trim(KEYWOR).Eq.'HFBTHOISON') Then Read(NFIREA,*) IF_THO,CBETHO counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IF_THO tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = CBETHO tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl End If ! Compute mass and charge of fission fragments If(trim(KEYWOR).Eq.'MASSFRAGME') Then Read(NFIREA,*) IFRAGM counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IFRAGM tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Compute and record densities of the fission fragments If(trim(KEYWOR).Eq.'FRAGDENSIT') Then Read(NFIREA,*) IDEALL,IDELOC,IDECON counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IDEALL tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IDELOC tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IDECON tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Rotate selected pairs of q.p. wave-functions If(trim(KEYWOR).Eq.'QPROTATION') Then Read(NFIREA,*) MIN_QP,DELTAE,XLOCMX,V2_MIN,ITRMAX,NTHETA counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = MIN_QP tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = DELTAE tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = XLOCMX tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = V2_MIN tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = ITRMAX tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = NTHETA tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If !=====================================! ! INTERACTIONS CONTROL FLAGS ! !=====================================! ! Switching to Energy Density Functional formalism If(trim(KEYWOR).Eq.'UNEDF_PROJ'.OR.trim(KEYWOR).Eq.'UNEDF-PROJ') Then Read(NFIREA,*) IF_EDF counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IF_EDF tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Skyrme force flags If(trim(KEYWOR).Eq.'SKYRME_STD'.OR.trim(KEYWOR).Eq.'SKYRME-STD') Then Read(NFIREA,*) ISTAND,KETA_J,KETA_W,KETACM,KETA_M counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = ISTAND tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = KETA_J tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = KETA_W tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = KETACM tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = KETA_M tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Type of Coulomb potential If(trim(KEYWOR).Eq.'COULOMBPAR') Then Read(NFIREA,*) ICOTYP,ICOUDI,ICOUEX counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = ICOTYP tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = ICOUDI tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = ICOUEX tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Gogny force If(trim(KEYWOR).Eq.'GOGNY') Then Read(NFIREA,*) I_GOGA counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = I_GOGA tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If !==============================! ! PAIRING CONTROL FLAGS ! !==============================! ! Gogny pairing If(trim(KEYWOR).Eq.'GOGNY_PAIR') Then Read(NFIREA,*) IGOGPA counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IGOGPA tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Switches on/off pairing If(trim(KEYWOR).Eq.'PAIRING') Then Read(NFIREA,*) IPAIRI counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IPAIRI tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! HFB pairing If(trim(KEYWOR).Eq.'HFB') Then Read(NFIREA,*) IPAHFB counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IPAHFB tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Defines the pairing window cut-off If(trim(KEYWOR).Eq.'CUTOFF') Then Read(NFIREA,*) ECUTOF counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = ECUTOF tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl End If ! Diagonalizes the mean-field (if HFB) If(trim(KEYWOR).Eq.'HFBMEANFLD') Then Read(NFIREA,*) IMFHFB counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IMFHFB tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Pairing interaction (neutrons only) If(trim(KEYWOR).Eq.'PAIRNINTER') Then Read(NFIREA,*) PRHO_N,PRHOSN,POWERN counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = PRHO_N tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = PRHOSN tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = POWERN tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl End If ! Pairing interaction (protons only) If(trim(KEYWOR).Eq.'PAIRPINTER') Then Read(NFIREA,*) PRHO_P,PRHOSP,POWERP counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = PRHO_P tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = PRHOSP tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = POWERP tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl End If ! Pairing interaction If(trim(KEYWOR).Eq.'PAIR_INTER') Then Read(NFIREA,*) PRHO_T,PRHOST,POWERT counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = PRHO_T tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = PRHOST tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = POWERT tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl End If ! Lipkin-Nogami correction If(trim(KEYWOR).Eq.'LIPKIN') Then Read(NFIREA,*) LIPKIN,LIPKIP counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = LIPKIN tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = LIPKIP tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Reading the fixed Lipkin-Nogami lambda2 parameters If(trim(KEYWOR).Eq.'FIXLAMB2_N'.OR.trim(KEYWOR).EQ.'FIXLAMB2-N') Then Read(NFIREA,*) FE2FIN,IF2FIN counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = FE2FIN tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IF2FIN tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If If(trim(KEYWOR).Eq.'FIXLAMB2_P'.OR.trim(KEYWOR).EQ.'FIXLAMB2-P') Then Read(NFIREA,*) FE2FIP,IF2FIP counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = FE2FIP tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IF2FIP tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If !==============================! ! SYMMETRY CONTROL FLAGS ! !==============================! ! time-reversal If(trim(KEYWOR).Eq.'ROTATION') Then Read(NFIREA,*) IROTAT counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IROTAT tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! y-simplex If(trim(KEYWOR).Eq.'SIMPLEXY') Then Read(NFIREA,*) ISIMPY counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = ISIMPY tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! y-signature If(trim(KEYWOR).Eq.'SIGNATUREY') Then Read(NFIREA,*) ISIGNY counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = ISIGNY tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! parity If(trim(KEYWOR).Eq.'PARITY') Then Read(NFIREA,*) IPARTY counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IPARTY tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! All 3 simplexes at once If(trim(KEYWOR).Eq.'TSIMPLEX3D') Then Read(NFIREA,*) ISIMTX,ISIMTY,ISIMTZ counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = ISIMTX tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = ISIMTY tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = ISIMTZ tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If !=====================================! ! HO BASIS CONTROL FLAGS ! !=====================================! ! Automatic adjustment of basis deformation and frequency If(trim(KEYWOR).Eq.'BASISAUTOM') Then Read(NFIREA,*) IBASIS counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IBASIS tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Flag controling the basis frequency If(trim(KEYWOR).Eq.'HOMEGAZERO') Then Read(NFIREA,*) FCHOM0 counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = FCHOM0 tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl End If ! Size of the HO basis If(trim(KEYWOR).Eq.'BASIS_SIZE'.OR.trim(KEYWOR).Eq.'BASIS-SIZE') Then Read(NFIREA,*) NOSCIL,NLIMIT,ENECUT If(NLIMIT.Eq.0) NLIMIT=((NOSCIL+1)*(NOSCIL+2)*(NOSCIL+3))/6 counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = NOSCIL tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = NLIMIT tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = ENECUT tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl End If ! Deformation of the HO basis If(trim(KEYWOR).Eq.'SURFAC_PAR'.OR.trim(KEYWOR).Eq.'SURFAC-PAR') Then Read(NFIREA,*) INNUMB,IZNUMB,R0PARM counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = INNUMB tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IZNUMB tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = R0PARM tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl End If ! Deformation of the HO basis If(trim(KEYWOR).Eq.'SURFAC_DEF'.OR.trim(KEYWOR).Eq.'SURFAC-DEF') Then counterDef = 0 LAMBDA=-1 Do While (LAMBDA < 0) Read(NFIREA,*) LAMBDA,MIU,ALPHAR counterDef = counterDef + 1 counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = ABS(LAMBDA) tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = MIU tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = ALPHAR tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl End Do counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = counterDef tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Automatic number of Gauss-Hermite points If(trim(KEYWOR).Eq.'OPTI_GAUSS'.OR.trim(KEYWOR).Eq.'OPTI-GAUSS') Then Read(NFIREA,*) IOPTGS counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IOPTGS tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Number of Gauss-Hermite points If(trim(KEYWOR).Eq.'GAUSHERMIT') Then Read(NFIREA,*) NXHERM,NYHERM,NZHERM counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = NXHERM tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = NYHERM tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = NZHERM tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If !=====================================! ! MULTIPOLE MOMENTS CONTROL FLAGS ! !=====================================! ! Linear constraints based on the RPA matrix If(trim(KEYWOR).Eq.'RPA_CONSTR'.OR.trim(KEYWOR).Eq.'RPA-CONSTR') Then Read(NFIREA,*) IF_RPA counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IF_RPA tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Flags controlling the number of constraints If(trim(KEYWOR).Eq.'MAX_MULTIP'.OR.trim(KEYWOR).Eq.'MAX-MULTIP') Then Read(NFIREA,*) NMUCON,NMUCOU,NMUPRI counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = NMUCON tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = NMUCOU tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = NMUPRI tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Deformations of the WS potential used to provide initial solution If(trim(KEYWOR).Eq.'DEFORMS_WS'.OR.trim(KEYWOR).Eq.'DEFORMS-WS') Then counterWS = 0 LAMBDA=-1 Do While (LAMBDA < 0) Read(NFIREA,*) LAMBDA,MIU counterWS = counterWS + 1 counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = ABS(LAMBDA) tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = MIU tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End Do counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = counterWS tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Additional deformations that will be used to optimize ! the WS deformation parameters (extra-push) If(trim(keywor).Eq.'PUSH_DEFOR'.OR.trim(keywor).Eq.'PUSH-DEFOR') Then L_PUSH = -1 N_PUSH = 0 Do While ( L_PUSH < 0 ) Read(NFIREA,*) L_PUSH, M_PUSH, Q_PUSH N_PUSH = N_PUSH + 1 counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = L_PUSH tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = M_PUSH tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = Q_PUSH tempNodeDbl%numero = counterDbl tempNodeDbl%name = keywor tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl End Do counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = N_PUSH tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Number of iterations to keep the extra constraints on If(trim(KEYWOR).Eq.'PUSH_ITERS'.OR.trim(KEYWOR).Eq.'PUSH-ITERS') Then Read(NFIREA,*) DOPUSH,NPLAST counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = DOPUSH tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = NPLAST tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Step of the minimization over WS deformations If(trim(KEYWOR).Eq.'STEP_MINIM'.OR.trim(KEYWOR).Eq.'STEP-MINIM') Then Read(NFIREA,*) STP_WS,ITE_WS counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = STP_WS tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = ITE_WS tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Choosing between a WS or HF-based mean-field If(trim(KEYWOR).Eq.'MEAN-FIELD'.OR.trim(KEYWOR).Eq.'MEAN_FIELD') Then Read(NFIREA,*) IHFPOT counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IHFPOT tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If !=====================================! ! ANGULAR MOMENTUM CONTROL FLAGS ! !=====================================! ! Rotational frequency If(trim(KEYWOR).Eq.'OMEGAY') Then Read(NFIREA,*) OMEGAY counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = OMEGAY tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl End If !================================! ! OUTPUT FILE CONTROL FLAGS ! !================================! ! Type of print-out If(trim(KEYWOR).Eq.'ONE_LINE'.OR.trim(KEYWOR).Eq.'ONE-LINE') Then Read(NFIREA,*) I1LINE counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = I1LINE tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Type of Nilsson label If(trim(KEYWOR).Eq.'NILSSONLAB') Then Read(NFIREA,*) NILXYZ counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = NILXYZ tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Flags controlling the deformations extracted from multipoles If(trim(KEYWOR).Eq.'BOHR_BETAS'.OR.trim(KEYWOR).Eq.'BOHR-BETAS') Then Read(NFIREA,*) NEXBET,IPRIBE,IPRIBL counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = NEXBET tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IPRIBE tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IPRIBL tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If !================================! ! I/O FLAGS CONTROL FLAGS ! !================================! ! Optionally recording ascii file If(trim(KEYWOR).Eq.'REVIEW') Then Read(NFIREA,*) IREVIE counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IREVIE tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Restart on densities in Gauss-Hermite mesh If(trim(KEYWOR).Eq.'FIELD_SAVE'.OR.trim(KEYWOR).Eq.'FIELD-SAVE') Then Read(NFIREA,*) IWRIFI counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IWRIFI tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If !==================================! ! CHECKPOINTING CONTROL FLAGS ! !==================================! ! Restart on densities in Gauss-Hermite mesh If(trim(KEYWOR).Eq.'RESTART') Then Read(NFIREA,*) ICONTI counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = ICONTI tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Restart on pairing If(trim(KEYWOR).Eq.'CONT_PAIRI'.OR.trim(KEYWOR).Eq.'CONT-PAIRI') Then Read(NFIREA,*) IPCONT counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IPCONT tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Restart on Lipkin-Nogami If(trim(KEYWOR).Eq.'CONTLIPKIN') Then Read(NFIREA,*) ILCONT counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = ILCONT tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Restart on HFB mean-fields If(trim(KEYWOR).Eq.'CONTFIELDS') Then Read(NFIREA,*) IFCONT counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IFCONT tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Restart on linear constraint If(trim(KEYWOR).Eq.'CONTAUGMEN') Then Read(NFIREA,*) IACONT counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IACONT tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Restart on Gogny mean-field If(trim(KEYWOR).Eq.'CONTGOGNY') Then Read(NFIREA,*) IGCONT counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IGCONT tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Restart on Gogny pairing-field If(trim(KEYWOR).Eq.'CONTGOGPAI') Then Read(NFIREA,*) IGPCON counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IGPCON tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Restart on WS binary file If(trim(KEYWOR).Eq.'READ_WOODS'.OR.trim(KEYWOR).Eq.'READ-WOODS') Then Read(NFIREA,*) IREAWS counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = IREAWS tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If !===================================! ! BACKWARD COMPATIBILTY WITH HFODD ! !===================================! If(trim(KEYWOR).Eq.'EXECUTE') Then NO_RUN=NO_RUN+1 End If End Do ! end of main DO WHILE loop Close(NFIREA) counterVector(1) = counterInt counterVector(2) = counterDbl counterVector(3) = counterStr End If ! myrank.eq.0 ! Broadcasting all the sizes (of integer, real and string data) Call MPI_Bcast(counterVector, 3, MPI_INTEGER, 0, MPI_COMM_WORLD, mpi_err) ! Reassigning local variables from broadcast values counterInt = counterVector(1) counterDbl = counterVector(2) counterStr = counterVector(3) Allocate(vectorInt_value(1:counterInt)) Allocate(vectorDbl_value(1:counterDbl)) Allocate(vectorStr_value(1:counterStr)) counterTot = counterInt + counterDbl + counterStr Allocate(vectorAll_names(1:counterTot)) ! Only processor 0 can fill out the vectors If(mpi_rank.Eq.0) Then ! Now that the final size of these vectors is known, we can ! allocate an array and broadcast to everyone else ! First, vector of integer data, and the relevant keywords Allocate(currentNodeInt) currentNodeInt => listeInt ii = 0 Do While ( associated(currentNodeInt) ) ii = ii + 1 vectorInt_value(ii) = currentNodeInt%value vectorAll_names(ii) = currentNodeInt%name Allocate(previousNodeInt) previousNodeInt => currentNodeInt currentNodeInt => currentNodeInt%next Deallocate(previousNodeInt) End Do ! Second, vector of real data, and the relevant keywords Allocate(currentNodeDbl) currentNodeDbl => listeDbl ii = 0 Do While ( associated(currentNodeDbl) ) ii = ii + 1 jj = ii + counterInt vectorDbl_value(ii) = currentNodeDbl%value vectorAll_names(jj) = currentNodeDbl%name Allocate(previousNodeDbl) previousNodeDbl => currentNodeDbl currentNodeDbl => currentNodeDbl%next Deallocate(previousNodeDbl) End Do ! Third: vector of string data, and the relevant keywords Allocate(currentNodeStr) currentNodeStr => listeStr ii = 0 Do While ( associated(currentNodeStr) ) ii = ii + 1 jj = ii + counterInt + counterDbl vectorStr_value(ii) = currentNodeStr%value vectorAll_names(jj) = currentNodeStr%name Allocate(previousNodeStr) previousNodeStr => currentNodeStr currentNodeStr => currentNodeStr%next Deallocate(previousNodeStr) End Do End If ! end of proc = 0 if ! Broadcasting all data of all types Call MPI_Bcast(vectorInt_value, counterInt, MPI_INTEGER, 0, MPI_COMM_WORLD, mpi_err) Call MPI_Bcast(vectorDbl_value, counterDbl, MPI_DOUBLE_PRECISION, 0, MPI_COMM_WORLD, mpi_err) Call MPI_Bcast(vectorStr_value, counterStr, mpi_string_10, 0, MPI_COMM_WORLD, mpi_err) ! Broadcasting all names (keywords needed to reconstruct data locally) Call MPI_Bcast(vectorAll_names, counterTot, mpi_string_10, 0, MPI_COMM_WORLD, mpi_err) ! Derived MPI type is not needed anymore, we release it Call mpi_type_free(mpi_string_10, mpi_err) End Subroutine mpi_getSequentialData !-------------------------------------------------------------------! ! This subroutine fills out the common blocks and variables that ! ! HFODD needs based on the data read in mpi_nameli() and broadcast ! ! by it. The subroutine is executed by all processors ! !-------------------------------------------------------------------! Subroutine mpi_setSequentialData(IOPTGS,NXHERM,NYHERM,NZHERM, & IPAIRI,IPAHFB,IPABCS,IMFHFB, & LIPKIN,LIPKIP, & FE2FIN,IF2FIN,FE2FIP,IF2FIP, & IROTAT,IREVIE, & ICONTI,IPCONT,IGCONT,IGPCON, & ILCONT,IFCONT,IACONT, & NUCHAO,NUPING,EPSPNG,EPSITE, & ISIMPY,ISIGNY,IPARTY, & ISIMTX,ISIMTY,ISIMTZ, & NOITER,NULAST, & NOSCIL,NLIMIT,ENECUT, & ICOTYP,ICOUDI,ICOUEX, & I_GOGA,IGOGPA, & I1LINE, & SLOWEV,SLOWOD,SLOWPA,SLOWLI, & ECUTOF, & NILXYZ, & IBROYD,N_ITER,BROTRI,ALPHAM, & PRHO_T,PRHOST,POWERT, & OMEGAY, & IFSHEL,MPOLYN,G_STUN,G_STUP,H_OMFA, & LAMBAS,MIUBAS,BETBAS,NUMDEF, & LAM_WS,MIU_WS,NUM_WS,ITE_WS,STP_WS, & LAMPUS,MIUPUS,NUMPUS,QMOPUS,IFPUSH, & NITPUS,IREAWS, & ISTAND,KETA_J,KETA_W,KETACM,KETA_M, & IF_RPA,IF_EDF,IF_THO,CBETHO, & IBASIS,IHFPOT,IFRAGM, & NMUCON,NMUCOU,NMUPRI, & NEXBET,IPRIBE,IPRIBL, & FCHOM0, & L_COLL,M_COLL,N_COLL, & IWRIFI, & INNUMB,IZNUMB,R0PARM, & DELTAE,XLOCMX,V2_MIN,ITRMAX,NTHETA, & MIN_QP,IDEALL,IDELOC,IDECON) Integer, intent(inout) :: ICOTYP, ICOUDI, ICOUEX, I_GOGA, IGOGPA, NOITER, IBROYD, & N_ITER, NULAST, NUPING, NUCHAO, ISIMPY, ISIGNY, IPARTY, & IROTAT, ISIMTX, ISIMTY, ISIMTZ, IPAIRI, IPAHFB, IPABCS, & IMFHFB, LIPKIN, LIPKIP, IOPTGS, NXHERM, NYHERM, NZHERM, & NOSCIL, NLIMIT, I1LINE, NILXYZ, IREVIE, ICONTI, IPCONT, & IGCONT, IGPCON, ILCONT, IFCONT, IACONT, IFSHEL, MPOLYN, & NUMDEF, NUM_WS, ITE_WS, IREAWS, ISTAND, KETA_J, KETA_W, & KETACM, KETA_M, IF_RPA, IF_EDF, IF_THO, IBASIS, IHFPOT, & IFRAGM, NUMPUS, IFPUSH, NITPUS, NMUCON, NMUCOU, NMUPRI, & NEXBET, IPRIBE, IPRIBL, N_COLL, IWRIFI, INNUMB, IZNUMB, & ITRMAX, NTHETA, MIN_QP, IDEALL, IDELOC, IDECON, IF2FIN, & IF2FIP Real, intent(inout) :: BROTRI, SLOWEV, SLOWOD, SLOWPA, SLOWLI, EPSITE, & EPSPNG, ECUTOF, ENECUT, PRHO_T, PRHOST, POWERT, & OMEGAY, G_STUN, G_STUP, H_OMFA, STP_WS, ALPHAM, & CBETHO, FCHOM0, R0PARM, DELTAE, XLOCMX, V2_MIN, & FE2FIN, FE2FIP Integer :: ii, jj, N_DUMM, counterDef, counterWS, counterColl, counter_push Integer, Dimension(:), pointer :: LAMBAS, MIUBAS, LAM_WS, MIU_WS, LAMPUS, MIUPUS Integer, Dimension(:), pointer :: L_COLL, M_COLL Real :: B_DUMM, PRHO_N, PRHOSN, PRHODN, POWERN, PRHO_P, PRHOSP, PRHODP, POWERP Real, Dimension(:), pointer :: BETBAS, QMOPUS Character(Len=4) :: SKYRME Character(Len=10) :: KEYWOR, string_10 ! common blocks for HFODD COMMON /CCPPAI/ PRHO_N,PRHODN,PRHOSN,POWERN, & PRHO_P,PRHODP,PRHOSP,POWERP ! Initialize local variables which are not initialized in PREDEF() MPOLYN=6; G_STUN=1.2; G_STUP=1.2; H_OMFA=3.5 NULLIFY(LAM_WS) NULLIFY(MIU_WS) NULLIFY(L_COLL) NULLIFY(M_COLL) NULLIFY(LAMPUS) NULLIFY(MIUPUS) NULLIFY(QMOPUS) NUM_WS=1 STP_WS=0.01 NULLIFY(LAMBAS) NULLIFY(MIUBAS) NULLIFY(BETBAS) NUMDEF=1 NUMPUS=0 !=================================================! ! ! ! RECONSTRUCTION OF THE INPUT OF TYPE INTEGER ! ! ! !=================================================! ii = 1 Do While ( ii .le. counterInt ) KEYWOR = vectorAll_names(ii) !==============================================! ! INTEGER: ITERATIONS CONTROL FLAGS ! !==============================================! ! Number of iterations If(trim(KEYWOR).Eq.'ITERATIONS') Then NOITER = vectorInt_value(ii) End If ! Switches on/off the Broyden method If(trim(KEYWOR).Eq.'BROYDEN') Then N_DUMM = vectorInt_value(ii) If(N_DUMM.GT.0) N_ITER=N_DUMM ii = ii +1 IBROYD = vectorInt_value(ii) End If ! Defines divergence If(trim(KEYWOR).Eq.'MAXANTIOSC') Then NULAST = vectorInt_value(ii) End If ! Control ping-pong effects If(trim(KEYWOR).Eq.'PING_PONG'.OR.trim(KEYWOR).Eq.'PING-PONG') Then NUPING = vectorInt_value(ii) End If ! Control chaotic divergence If(trim(KEYWOR).Eq.'CHAOTIC') Then NUCHAO = vectorInt_value(ii) End If !==============================================! ! INTEGER: SPECIFIC FEATURES CONTROL FLAGS ! !==============================================! ! Type of shell correction If(trim(KEYWOR).Eq.'SHELLCORCT') Then IFSHEL = vectorInt_value(ii) End If ! Number of polynomials for the shell correction If(trim(KEYWOR).Eq.'SHELLPARAM') Then MPOLYN = vectorInt_value(ii) End If ! Switching to interface with HFBTHO If(trim(KEYWOR).Eq.'HFBTHOISON') Then IF_THO = vectorInt_value(ii) End If ! Compute mass and charge of fission fragments If(trim(KEYWOR).Eq.'MASSFRAGME') Then IFRAGM = vectorInt_value(ii) End If ! Compute and record densities of the fission fragments If(trim(KEYWOR).Eq.'FRAGDENSIT') Then IDECON = vectorInt_value(ii) ii = ii + 1 IDELOC = vectorInt_value(ii) ii = ii + 1 IDEALL = vectorInt_value(ii) End If ! Rotate selected pairs of q.p. wave-functions If(trim(KEYWOR).Eq.'QPROTATION') Then NTHETA = vectorInt_value(ii) ii = ii + 1 ITRMAX = vectorInt_value(ii) ii = ii + 1 MIN_QP = vectorInt_value(ii) End If ! Definition of the collective space used for calculation of the ! collective inertia If(trim(KEYWOR).Eq.'COLL_SPACE'.OR.trim(KEYWOR).Eq.'COLL-SPACE') Then counterColl = vectorInt_value(ii) N_COLL = counterColl Allocate(L_COLL(1:counterColl)) Allocate(M_COLL(1:counterColl)) Do jj = 1, counterColl ii = ii + 1 M_COLL(jj) = vectorInt_value(ii) ii = ii + 1 L_COLL(jj) = vectorInt_value(ii) End Do End If !==============================================! ! INTEGER: INTERACTIONS CONTROL FLAGS ! !==============================================! ! Switching to Energy Density Functional formalism If(trim(KEYWOR).Eq.'UNEDF_PROJ'.OR.trim(KEYWOR).Eq.'UNEDF-PROJ') Then IF_EDF = vectorInt_value(ii) End If ! Skyrme force flags If(trim(KEYWOR).Eq.'SKYRME_STD'.OR.trim(KEYWOR).Eq.'SKYRME-STD') Then KETA_M = vectorInt_value(ii) ii = ii + 1 KETACM = vectorInt_value(ii) ii = ii + 1 KETA_W = vectorInt_value(ii) ii = ii + 1 KETA_J = vectorInt_value(ii) ii = ii + 1 ISTAND = vectorInt_value(ii) End If ! Type of Coulomb potential If(trim(KEYWOR).Eq.'COULOMBPAR') Then ICOUEX = vectorInt_value(ii) ii = ii + 1 ICOUDI = vectorInt_value(ii) ii = ii + 1 ICOTYP = vectorInt_value(ii) End If ! Gogny force If(trim(KEYWOR).Eq.'GOGNY') Then I_GOGA = vectorInt_value(ii) End If !==============================================! ! INTEGER: PAIRING CONTROL FLAGS ! !==============================================! ! Gogny pairing If(trim(KEYWOR).Eq.'GOGNY_PAIR'.OR.trim(KEYWOR).Eq.'GOGNY-PAIR') Then IGOGPA = vectorInt_value(ii) End If ! Switches on/off pairing If(trim(KEYWOR).Eq.'PAIRING') Then IPAIRI = vectorInt_value(ii) End If ! HFB pairing If(trim(KEYWOR).Eq.'HFB') Then IPAHFB = vectorInt_value(ii) End If ! BCS pairing If(trim(KEYWOR).Eq.'BCS') Then IPABCS = vectorInt_value(ii) End If ! Diagonalizes the mean-field (if HFB) If(trim(KEYWOR).Eq.'HFBMEANFLD') Then IMFHFB = vectorInt_value(ii) End If ! Lipkin-Nogami correction If(trim(KEYWOR).Eq.'LIPKIN') Then LIPKIP = vectorInt_value(ii) ii = ii + 1 LIPKIN = vectorInt_value(ii) End If ! Reading the fixed Lipkin-Nogami lambda2 parameters If(trim(KEYWOR).Eq.'FIXLAMB2_N'.OR.trim(KEYWOR).EQ.'FIXLAMB2-N') Then IF2FIN = vectorInt_value(ii) End If If(trim(KEYWOR).Eq.'FIXLAMB2_P'.OR.trim(KEYWOR).EQ.'FIXLAMB2-P') Then IF2FIP = vectorInt_value(ii) End If !==============================================! ! INTEGER: SYMMETRY CONTROL FLAGS ! !==============================================! ! time-reversal If(trim(KEYWOR).Eq.'ROTATION') Then IROTAT = vectorInt_value(ii) End If ! y-simplex If(trim(KEYWOR).Eq.'SIMPLEXY') Then ISIMPY = vectorInt_value(ii) End If ! y-signature If(trim(KEYWOR).Eq.'SIGNATUREY') Then ISIGNY = vectorInt_value(ii) End If ! parity If(trim(KEYWOR).Eq.'PARITY') Then IPARTY = vectorInt_value(ii) End If ! All 3 simplexes at once If(trim(KEYWOR).Eq.'TSIMPLEX3D') Then ISIMTZ = vectorInt_value(ii) ii = ii + 1 ISIMTY = vectorInt_value(ii) ii = ii + 1 ISIMTX = vectorInt_value(ii) End If !==============================================! ! INTEGER: HO BASIS CONTROL FLAGS ! !==============================================! ! Automatic adjustment of HO basis parameters If(trim(KEYWOR).Eq.'BASISAUTOM') Then IBASIS = vectorInt_value(ii) End If ! Size of the HO basis If(trim(KEYWOR).Eq.'BASIS_SIZE'.OR.trim(KEYWOR).Eq.'BASIS-SIZE') Then NLIMIT = vectorInt_value(ii) ii = ii + 1 NOSCIL = vectorInt_value(ii) End If ! Particle numbers of the HO basis If(trim(KEYWOR).Eq.'SURFAC_PAR'.OR.trim(KEYWOR).Eq.'SURFAC-PAR') Then IZNUMB = vectorInt_value(ii) ii = ii + 1 INNUMB = vectorInt_value(ii) End If ! Deformation of the HO basis If(trim(KEYWOR).Eq.'SURFAC_DEF'.OR.trim(KEYWOR).Eq.'SURFAC-DEF') Then counterDef = vectorInt_value(ii) NUMDEF = counterDef Allocate(LAMBAS(1:counterDef)) Allocate(MIUBAS(1:counterDef)) Do jj = 1, counterDef ii = ii + 1 MIUBAS(jj) = vectorInt_value(ii) ii = ii + 1 LAMBAS(jj) = vectorInt_value(ii) End Do End If ! Automatic number of Gauss-Hermite points If(trim(KEYWOR).Eq.'OPTI_GAUSS'.OR.trim(KEYWOR).Eq.'OPTI-GAUSS') Then IOPTGS = vectorInt_value(ii) End If ! Number of Gauss-Hermite points If(trim(KEYWOR).Eq.'GAUSHERMIT') Then NZHERM = vectorInt_value(ii) ii = ii + 1 NYHERM = vectorInt_value(ii) ii = ii + 1 NXHERM = vectorInt_value(ii) End If !==============================================! ! INTEGER: MULTIPOLE MOMENT CONTROL FLAGS ! !==============================================! If(trim(KEYWOR).Eq.'RPA_CONSTR'.OR.trim(KEYWOR).Eq.'RPA-CONSTR') Then IF_RPA = vectorInt_value(ii) End If ! Number of iterations to keep the extra constraints on If(trim(KEYWOR).Eq.'MAX_MULTIP'.OR.trim(KEYWOR).Eq.'MAX-MULTIP') Then NMUPRI = vectorInt_value(ii) ii = ii + 1 NMUCOU = vectorInt_value(ii) ii = ii + 1 NMUCON = vectorInt_value(ii) End If ! Deformations of the WS potential used to provide initial solution If(trim(KEYWOR).Eq.'DEFORMS_WS'.OR.trim(KEYWOR).Eq.'DEFORMS-WS') Then counterWS = vectorInt_value(ii) NUM_WS = counterWS Allocate(LAM_WS(1:counterWS)) Allocate(MIU_WS(1:counterWS)) Do jj = 1, counterWS ii = ii + 1 MIU_WS(jj) = vectorInt_value(ii) ii = ii + 1 LAM_WS(jj) = vectorInt_value(ii) End Do End If ! Additional deformations that will be used to optimize ! the WS deformation parameters (extra-push) If(trim(KEYWOR).Eq.'PUSH_DEFOR'.OR.trim(KEYWOR).Eq.'PUSH-DEFOR') Then counter_push = vectorInt_value(ii) NUMPUS = counter_push Allocate(LAMPUS(1:counter_push)) Allocate(MIUPUS(1:counter_push)) Do jj = 1, counter_push ii = ii + 1 MIUPUS(jj) = vectorInt_value(ii) ii = ii + 1 LAMPUS(jj) = vectorInt_value(ii) End Do End If ! Number of iterations to keep the extra constraints on If(trim(KEYWOR).Eq.'PUSH_ITERS'.OR.trim(KEYWOR).Eq.'PUSH-ITERS') Then NITPUS = vectorInt_value(ii) ii = ii + 1 IFPUSH = vectorInt_value(ii) End If ! Number of iterations of the minimization If(trim(KEYWOR).Eq.'STEP_MINIM'.OR.trim(KEYWOR).Eq.'STEP-MINIM') Then ITE_WS = vectorInt_value(ii) End If ! Choosing between a WS or HF-based mean-field If(trim(KEYWOR).Eq.'MEAN-FIELD'.OR.trim(KEYWOR).Eq.'MEAN_FIELD') Then IHFPOT = vectorInt_value(ii) End If !==============================================! ! INTEGER: OUTPUT FILE CONTROL FLAGS ! !==============================================! ! Type of print-out If(trim(KEYWOR).Eq.'ONE_LINE'.OR.trim(KEYWOR).Eq.'ONE-LINE') Then I1LINE = vectorInt_value(ii) End If ! Type of Nilsson label If(trim(KEYWOR).Eq.'NILSSONLAB') Then NILXYZ = vectorInt_value(ii) End If ! Number of iterations to keep the extra constraints on If(trim(KEYWOR).Eq.'BOHR_BETAS'.OR.trim(KEYWOR).Eq.'BOHR-BETAS') Then IPRIBL = vectorInt_value(ii) ii = ii + 1 IPRIBE = vectorInt_value(ii) ii = ii + 1 NEXBET = vectorInt_value(ii) End If !==============================================! ! INTEGER: I/O FLAGS CONTROL FLAGS ! !==============================================! ! Optionally recording ascii file If(trim(KEYWOR).Eq.'REVIEW') Then IREVIE = vectorInt_value(ii) End If ! Restart on densities in Gauss-Hermite mesh If(trim(KEYWOR).Eq.'FIELD_SAVE'.OR.trim(KEYWOR).Eq.'FIELD-SAVE') Then IWRIFI = vectorInt_value(ii) End If !==============================================! ! INTEGER: CHECKPOINTING CONTROL FLAGS ! !==============================================! ! Restart on densities If(trim(KEYWOR).Eq.'RESTART') Then ICONTI = vectorInt_value(ii) End If ! Restart on pairing If(trim(KEYWOR).Eq.'CONT_PAIRI'.OR.trim(KEYWOR).Eq.'CONT-PAIRI') Then IPCONT = vectorInt_value(ii) End If ! Restart on Lipkin-Nogami If(trim(KEYWOR).Eq.'CONTLIPKIN') Then ILCONT = vectorInt_value(ii) End If ! Restart on HFB mean-fields If(trim(KEYWOR).Eq.'CONTFIELDS') Then IFCONT = vectorInt_value(ii) End If ! Restart on linear constraint If(trim(KEYWOR).Eq.'CONTAUGMEN') Then IACONT = vectorInt_value(ii) End If ! Restart on Gogny mean-field If(trim(KEYWOR).Eq.'CONTGOGNY') Then IGCONT = vectorInt_value(ii) End If ! Restart on Gogny pairing-field If(trim(KEYWOR).Eq.'CONTGOGPAI') Then IGPCON = vectorInt_value(ii) End If ! Restart on WS binary files If(trim(KEYWOR).Eq.'READ_WOODS'.OR.trim(KEYWOR).Eq.'READ-WOODS') Then IREAWS = vectorInt_value(ii) End If ! Go to the next item on the list ii = ii + 1 End Do Deallocate(vectorInt_value) !==============================================! ! ! ! RECONSTRUCTION OF THE INPUT OF TYPE REAL ! ! ! !==============================================! ii = 1 Do While ( ii .le. counterDbl ) KEYWOR = vectorAll_names(ii + counterInt) !===========================================! ! REAL: ITERATIONS CONTROL FLAGS ! !===========================================! ! Switches on/off the Broyden method If(trim(KEYWOR).Eq.'BROYDEN') Then B_DUMM = vectorDbl_value(ii) If(B_DUMM.GT.0.0D0) BROTRI=B_DUMM ii = ii + 1 ALPHAM = vectorDbl_value(ii) End If ! Slow-down parameters for the mean-field If(trim(KEYWOR).Eq.'SLOW_DOWN'.OR.trim(KEYWOR).Eq.'SLOW-DOWN') Then SLOWOD = vectorDbl_value(ii) ii = ii + 1 SLOWEV = vectorDbl_value(ii) End If ! Slow-down parameters for the pairing field If(trim(KEYWOR).Eq.'SLOW_PAIR'.OR.trim(KEYWOR).Eq.'SLOW-PAIR') Then SLOWPA = vectorDbl_value(ii) End If ! Slow-down parameters for the Lipkin-Nogami correction If(trim(KEYWOR).Eq.'SLOWLIPKIN') Then SLOWLI = vectorDbl_value(ii) End If ! Defines the precision of convergence If(trim(KEYWOR).Eq.'ITERAT_EPS'.OR.trim(KEYWOR).Eq.'ITERAT-EPS') Then EPSITE = vectorDbl_value(ii) End If ! Control ping-pong effects If(trim(KEYWOR).Eq.'PING_PONG'.OR.trim(KEYWOR).Eq.'PING-PONG') Then EPSPNG = vectorDbl_value(ii) End If !===========================================! ! REAL: SPECIFIC FEATURES CONTROL FLAGS ! !===========================================! ! Parameters for the shell correction If(trim(KEYWOR).Eq.'SHELLPARAM') Then H_OMFA = vectorDbl_value(ii) ii = ii + 1 G_STUP = vectorDbl_value(ii) ii = ii + 1 G_STUN = vectorDbl_value(ii) End If ! Switching to interface with HFBTHO If(trim(KEYWOR).Eq.'HFBTHOISON') Then CBETHO = vectorDbl_value(ii) End If ! Rotate selected pairs of q.p. wave-functions If(trim(KEYWOR).Eq.'QPROTATION') Then V2_MIN = vectorDbl_value(ii) ii = ii + 1 XLOCMX = vectorDbl_value(ii) ii = ii + 1 DELTAE = vectorDbl_value(ii) End If !===========================================! ! REAL: PAIRING CONTROL FLAGS ! !===========================================! ! Control chaotic divergence If(trim(KEYWOR).Eq.'CUTOFF') Then ECUTOF = vectorDbl_value(ii) End If ! Pairing Interaction If(trim(KEYWOR).Eq.'PAIR_INTER'.OR.trim(KEYWOR).Eq.'PAIR-INTER') Then POWERT = vectorDbl_value(ii) ii = ii + 1 PRHOST = vectorDbl_value(ii) ii = ii + 1 PRHO_T = vectorDbl_value(ii) End If ! The fixed Lipkin-Nogami lambda2 parameters If(trim(KEYWOR).Eq.'FIXLAMB2_N'.OR.trim(KEYWOR).Eq.'FIXLAMB2-N') Then FE2FIN = vectorDbl_value(ii) End If If(trim(KEYWOR).Eq.'FIXLAMB2_P'.OR.trim(KEYWOR).Eq.'FIXLAMB2-P') Then FE2FIP = vectorDbl_value(ii) End If ! Pairing Interaction (neutrons only) If(trim(KEYWOR).Eq.'PAIRNINTER') Then POWERN = vectorDbl_value(ii) ii = ii + 1 PRHOSN = vectorDbl_value(ii) ii = ii + 1 PRHO_N = vectorDbl_value(ii) End If ! Pairing Interaction (protons only) If(trim(KEYWOR).Eq.'PAIRPINTER') Then POWERP = vectorDbl_value(ii) ii = ii + 1 PRHOSP = vectorDbl_value(ii) ii = ii + 1 PRHO_P = vectorDbl_value(ii) End If !===========================================! ! REAL: HO BASIS CONTROL FLAGS ! !===========================================! ! Flag controling the basis frequency If(trim(KEYWOR).Eq.'HOMEGAZERO') Then FCHOM0 = vectorDbl_value(ii) End If ! Size of the HO basis If(trim(KEYWOR).Eq.'BASIS_SIZE'.OR.trim(KEYWOR).Eq.'BASIS-SIZE') Then ENECUT = vectorDbl_value(ii) End If ! Particle numbers of the HO basis If(trim(KEYWOR).Eq.'SURFAC_PAR'.OR.trim(KEYWOR).Eq.'SURFAC-PAR') Then R0PARM = vectorDbl_value(ii) End If ! Deformation of the HO basis If(trim(KEYWOR).Eq.'SURFAC_DEF'.OR.trim(KEYWOR).Eq.'SURFAC-DEF') Then Allocate(BETBAS(1:counterDef)) Do jj = 1, counterDef BETBAS(jj) = vectorDbl_value(ii) ii = ii + 1 End Do End If !===========================================! ! REAL: MULTIPOLE MOMENT CONTROL FLAGS ! !===========================================! ! Deformations of the WS potential If(trim(KEYWOR).Eq.'PUSH_DEFOR'.OR.trim(KEYWOR).Eq.'PUSH-DEFOR') Then Allocate(QMOPUS(1:counter_push)) DO jj = 1, counter_push QMOPUS(jj) = vectorDbl_value(ii) ii = ii + 1 End Do ii = ii - 1 End If ! Step of the minimization over WS deformations If(trim(KEYWOR).Eq.'STEP_MINIM'.OR.trim(KEYWOR).Eq.'STEP-MINIM') Then STP_WS = vectorDbl_value(ii) End If !===========================================! ! REAL: ANGULAR MOMENTUM CONTROL FLAGS ! !===========================================! ! Rotational frequency If(trim(KEYWOR).Eq.'OMEGAY') Then OMEGAY = vectorDbl_value(ii) End If ! Go to the next item on the list ii = ii + 1 End Do Deallocate(vectorDbl_value) !================================================! ! ! ! RECONSTRUCTION OF THE INPUT OF TYPE STRING*10 ! ! ! !================================================! ii = 1 Do While ( ii .le. counterStr ) KEYWOR = vectorAll_names(ii + counterInt + counterDbl) If(trim(KEYWOR).Eq.'SKYRME-SET') Then string_10 = vectorStr_value(ii) SKYRME = string_10(1:4) End If ! Go to the next item on the list ii = ii + 1 End Do Deallocate(vectorStr_value) Deallocate(vectorAll_names) End Subroutine mpi_setSequentialData !---------------------------------------------------------------! ! This SUBROUTINE reads data from input file. These data will ! ! then be used to construct a collection of jobs for HFODD. ! !---------------------------------------------------------------! Subroutine mpi_getParallelData() Integer, Parameter :: ndatin = 83 Character(Len=4) :: SKYRME Character(Len=80) :: filein Character(Len=10) :: keywor Integer :: ii,jj,ierror,numero,mpi_rank,mpi_err,mpi_string_10,mpidef,mpinck,mpibas Integer :: lambda,miu,numberQ,iq,lambda_res,miu_res,numberQ_res,iq_res,iq_batch Integer :: izstrt,izstep,nstepz,instrt,instep,nstepn,n_mini,n_step,noffre,nofshl, & nofb20,nsmini,nsstep,nofsta,ifcons,ibatch,nbatch,if_opt,iftemp,num_qn, & num_qn_res Real :: xzstrt,xzstep,xnstrt,xnstep,o_mini,o_step,b20min,b20stp, & qBegin,qFin,qBegin_res,qFin_res,qn_min,qn_max,qn_min_res, & qn_max_res,TEMP_T Type(Node_int), pointer :: tempNodeInt, currentNodeInt, previousNodeInt Type(Node_dbl), pointer :: tempNodeDbl, currentNodeDbl, previousNodeDbl Type(Node_str), pointer :: tempNodeStr, currentNodeStr, previousNodeStr ! Getting the rank of the current process Call mpi_comm_rank(MPI_COMM_WORLD, mpi_rank, mpi_err) ! Define MPI type for a character string of length stringLength and commit Call mpi_type_contiguous(stringLength, MPI_CHARACTER, mpi_string_10, mpi_err) Call mpi_type_commit(mpi_string_10, mpi_err) ! Only the root process reads the data If(mpi_rank.Eq.0) Then filein = "hfodd_mpiio.d" Open(ndatin,file=filein,status='old',form='formatted',iostat=ierror) If(ierror .NE. 0) Then Write(6,'("error in opening file : ",a80)') filein write(6,'("error code : ",i12)') ierror Stop 'error in readin - i/o' End If keywor="Initialize" counterInt = 0 counterDbl = 0 counterStr = 0 counterInt = counterInt + 1 Allocate(listeInt) listeInt%value = 0 listeInt%numero = counterInt listeInt%name = keywor NULLIFY(listeInt%next) counterDbl = counterDbl + 1 Allocate(listeDbl) listeDbl%value = 0.0 listeDbl%numero = counterDbl listeDbl%name = keywor NULLIFY(listeDbl%next) counterStr = counterStr + 1 Allocate(listeStr) listeStr%value = "INITIALIZE" listeStr%numero = counterStr listeStr%name = KEYWOR NULLIFY(listeStr%next) Do While (trim(keywor) .NE. 'END_DATA') Read(ndatin, FMT="(A)") keywor ! Type of calculation (deformation, deformation with restart, ! basis) If(trim(keywor) == 'CALCULMODE') Then Read(ndatin,*) mpidef, mpibas counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = mpidef tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = mpibas tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Setting up constrained or unconstrained calculations If(trim(keywor) == 'CONSTR_DEF'.OR.trim(keywor) == 'CONSTR-DEF') Then Read(ndatin,*) ifcons counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = ifcons tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Setting up batch mode for constrained calculations If(trim(keywor) == 'BATCH_MODE'.OR.trim(keywor) == 'BATCH-MODE') Then Read(ndatin,*) ibatch, nbatch counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = ibatch tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = nbatch tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Specifying which multipoles will be 'batched' If(trim(keywor) == 'BATCH_SPEC'.OR.trim(keywor) == 'BATCH-SPEC') Then lambda = -1 iq_batch = 0 Do While ( lambda < 0 ) iq_batch = iq_batch + 1 Read(ndatin,*) lambda, miu counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = lambda tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = miu tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt End Do counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = iq_batch tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Defining nuclear grid If(trim(keywor) == 'ALL_NUCLEI'.OR.trim(keywor) == 'ALL-NUCLEI') Then Read(ndatin,*) izstrt, izstep, nstepz, instrt, instep, nstepn counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = izstrt tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = izstep tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = nstepz tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = instrt tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = instep tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = nstepn tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Defining nuclear grid If(trim(keywor) == 'REALNUCLEI') Then Read(ndatin,*) xzstrt, xzstep, nstepz, xnstrt, xnstep, nstepn counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = xzstrt tempNodeDbl%numero = counterDbl tempNodeDbl%name = keywor tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = xzstep tempNodeDbl%numero = counterDbl tempNodeDbl%name = keywor tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = nstepz tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = xnstrt tempNodeDbl%numero = counterDbl tempNodeDbl%name = keywor tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = xnstep tempNodeDbl%numero = counterDbl tempNodeDbl%name = keywor tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = nstepn tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Defining default multipole grid If(trim(keywor) == 'MULTICONST') Then lambda = -1 iq = 0 Do While ( lambda < 0 ) iq = iq + 1 Read(ndatin,*) lambda, miu, qBegin, qFin, numberQ counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = abs(lambda) tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = miu tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = numberQ tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = qBegin tempNodeDbl%numero = counterDbl tempNodeDbl%name = keywor tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = qFin tempNodeDbl%numero = counterDbl tempNodeDbl%name = keywor tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl End Do counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = iq tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Defining new multipole grid (active only if mpidef=2). If ! active, information under keyword MULTICONST will be used ! to identify the restart points If(trim(keywor) == 'MULTIRESTA') Then lambda_res = -1 iq_res = 0 Do While ( lambda_res < 0 ) iq_res = iq_res + 1 Read(ndatin,*) lambda_res, miu_res, qBegin_res, qFin_res, numberQ_res counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = abs(lambda_res) tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = miu_res tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = numberQ_res tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = qBegin_res tempNodeDbl%numero = counterDbl tempNodeDbl%name = keywor tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = qFin_res tempNodeDbl%numero = counterDbl tempNodeDbl%name = keywor tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl End Do counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = iq_res tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Sets constraints on the size of the neck If(trim(keywor) == 'CALCULNECK') Then Read(ndatin,*) mpinck counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = mpinck tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Defining default grid of neck constraints If(trim(keywor) == 'NECK_CONST') Then Read(ndatin,*) qn_min, qn_max, num_qn counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = qn_min tempNodeDbl%numero = counterDbl tempNodeDbl%name = keywor tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = qn_max tempNodeDbl%numero = counterDbl tempNodeDbl%name = keywor tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = num_qn tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Defining new neck grid (active only if mpidef=2). If ! active, information under keyword NECK_CONST will be used ! to identify the restart points If(trim(keywor) == 'NECK_RESTA') Then Read(ndatin,*) qn_min_res, qn_max_res, num_qn_res counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = qn_min_res tempNodeDbl%numero = counterDbl tempNodeDbl%name = keywor tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = qn_max_res tempNodeDbl%numero = counterDbl tempNodeDbl%name = keywor tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = num_qn_res tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Choose to optimize the deformation grid by centering new ! multipolarity constraints on the restart equilibrium value If(trim(keywor) == 'OPTIM_GRID'.OR.trim(keywor) == 'OPTIM-GRID') Then Read(ndatin,*) if_opt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = if_opt tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Skyrme force name If(trim(KEYWOR).Eq.'ALL_FORCES'.OR.trim(KEYWOR).Eq.'ALL-FORCES') Then counterForce = 0 numero = -1 Do While (numero < 0) Read(ndatin,*) numero, SKYRME counterForce = counterForce + 1 counterStr = counterStr + 1 Allocate(tempNodeStr) tempNodeStr%value = trim(SKYRME)//" " tempNodeStr%numero = counterStr tempNodeStr%name = KEYWOR tempNodeStr%next => listeStr listeStr => tempNodeStr End Do counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = counterForce tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Deciding if the restart will operate on a grid of temperatures If(trim(keywor) == 'TEMPRESTRT') Then Read(ndatin,*) iftemp counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = iftemp tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! Choosing a grid of temperatures (mpidef>0 only) If(trim(keywor) == 'TEMPVALUES') Then counterTemperature = 0 numero = -1 Do While (numero < 0) Read(ndatin,*) numero, TEMP_T counterTemperature = counterTemperature + 1 counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = TEMP_T tempNodeDbl%numero = counterDbl tempNodeDbl%name = KEYWOR tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl End Do counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = counterTemperature tempNodeInt%numero = counterInt tempNodeInt%name = KEYWOR tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! For mpibas=1, grid for the basis frequency If(trim(keywor) == 'BASIS_FREQ'.OR.trim(keywor) == 'BASIS-FREQ') Then Read(ndatin,*) o_mini, o_step, noffre counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = o_mini tempNodeDbl%numero = counterDbl tempNodeDbl%name = keywor tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = o_step tempNodeDbl%numero = counterDbl tempNodeDbl%name = keywor tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = noffre tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! For mpibas=1, grid for the basis number of shells If(trim(keywor) == 'BASIS_NSHL'.OR.trim(keywor) == 'BASIS-NSHL') Then Read(ndatin,*) n_mini, n_step, nofshl counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = n_mini tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = n_step tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = nofshl tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! For mpibas=1, grid for the basis number of states If(trim(keywor) == 'BASIS_NSTA'.OR.trim(keywor) == 'BASIS-NSTA') Then Read(ndatin,*) nsmini, nsstep, nofsta counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = nsmini tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = nsstep tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = nofsta tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt End If ! For mpibas=1, grid for the basis deformation If(trim(keywor) == 'BASIS_DEFS'.OR.trim(keywor) == 'BASIS-DEFS') Then Read(ndatin,*) b20min, b20stp, nofb20 counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = b20min tempNodeDbl%numero = counterDbl tempNodeDbl%name = keywor tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterDbl = counterDbl + 1 Allocate(tempNodeDbl) tempNodeDbl%value = b20stp tempNodeDbl%numero = counterDbl tempNodeDbl%name = keywor tempNodeDbl%next => listeDbl listeDbl => tempNodeDbl counterInt = counterInt + 1 Allocate(tempNodeInt) tempNodeInt%value = nofb20 tempNodeInt%numero = counterInt tempNodeInt%name = keywor tempNodeInt%next => listeInt listeInt => tempNodeInt End If If(trim(keywor) == 'END_DATA') Then Exit End If End Do Close(ndatin) counterVector(1) = counterInt counterVector(2) = counterDbl counterVector(3) = counterStr End If ! Broadcasting all sizes of the different types of data Call MPI_Bcast(counterVector, 3, MPI_INTEGER, 0, MPI_COMM_WORLD, mpi_err) ! Reconstructing the counters locally and allocating local arrays counterInt = counterVector(1) counterDbl = counterVector(2) counterStr = counterVector(3) Allocate(vectorInt_value(1:counterInt)) Allocate(vectorDbl_value(1:counterDbl)) Allocate(vectorStr_value(1:counterStr)) ! Allocating a local array for data names counterTot = counterInt + counterDbl + counterStr Allocate(vectorAll_names(1:counterTot)) ! Only processor 0 can fill out the vectors If(mpi_rank.Eq.0) Then ! Now that the final size of these vectors is known, we can ! allocate an array and broadcast to everyone else ! First, vector of integer data, and the relevant keywords Allocate(currentNodeInt) currentNodeInt => listeInt ii = 0 Do While ( associated(currentNodeInt) ) ii = ii + 1 vectorInt_value(ii) = currentNodeInt%value vectorAll_names(ii) = currentNodeInt%name Allocate(previousNodeInt) previousNodeInt => currentNodeInt currentNodeInt => currentNodeInt%next Deallocate(previousNodeInt) End Do ! Second, vector of real data, and the relevant keywords Allocate(currentNodeDbl) currentNodeDbl => listeDbl ii = 0 Do While ( associated(currentNodeDbl) ) ii = ii + 1 jj = ii + counterInt vectorDbl_value(ii) = currentNodeDbl%value vectorAll_names(jj) = currentNodeDbl%name Allocate(previousNodeDbl) previousNodeDbl => currentNodeDbl currentNodeDbl => currentNodeDbl%next Deallocate(previousNodeDbl) End Do ! Third: vector of string data, and the relevant keywords Allocate(currentNodeStr) currentNodeStr => listeStr ii = 0 Do While ( associated(currentNodeStr) ) ii = ii + 1 jj = ii + counterInt + counterDbl vectorStr_value(ii) = currentNodeStr%value vectorAll_names(jj) = currentNodeStr%name Allocate(previousNodeStr) previousNodeStr => currentNodeStr currentNodeStr => currentNodeStr%next Deallocate(previousNodeStr) End Do End If ! end of proc = 0 if ! Broadcasting all data Call MPI_Bcast(vectorInt_value, counterInt, MPI_INTEGER, 0, MPI_COMM_WORLD, mpi_err) Call MPI_Bcast(vectorDbl_value, counterDbl, MPI_DOUBLE_PRECISION, 0, MPI_COMM_WORLD, mpi_err) Call MPI_Bcast(vectorStr_value, counterStr, mpi_string_10, 0, MPI_COMM_WORLD, mpi_err) ! Broadcasting all data names Call MPI_Bcast(vectorAll_names, counterTot, mpi_string_10, 0, MPI_COMM_WORLD, mpi_err) ! Releasing derived MPI datatype Call mpi_type_free(mpi_string_10, mpi_err) End Subroutine mpi_getParallelData !-------------------------------------------------------------------! ! This subroutine initializes processor-dependent variables used to ! ! construct the list of jobs. The subroutine is executed by all ! ! processors ! !-------------------------------------------------------------------! Subroutine mpi_initializeParallelData() Integer :: iq modeMPI_deformation = 1 modeMPI_neck = 0 modeMPI_basis = 0 turnConstraintOn = 1 optimizePath = 0 batch_mode = 0 number_batch = 1 number_constraint_batch = 1 counterForce = 1 Allocate(forceVector(1:counterForce)) forceVector(:) = 'SKM*' Allocate(lambda_batch(1:number_constraint_batch)) Allocate( miu_batch(1:number_constraint_batch)) lambda_batch(:) = 3 miu_batch(:) = 0 proton_initValue = 66; proton_stepValue = 2; proton_stepNumber = 1 neutron_initValue = 86; neutron_stepValue = 2; neutron_stepNumber = 2 x_proton_initValue = -1.0D0; x_proton_stepValue = 0.1 x_neutron_initValue = -1.0D0; x_neutron_stepValue = 0.1 restartTemperature = 0 counterTemperature = 1 Allocate(temperatureVector(1:counterTemperature)) temperatureVector(:) = 0.0D0 numberConstraints = 1 Allocate(qStepNumber(1:numberConstraints)) Allocate( qLambda(1:numberConstraints)) Allocate( qMiu(1:numberConstraints)) Allocate( qinitValue(1:numberConstraints)) Allocate( qFinValue(1:numberConstraints)) Do iq = 1, numberConstraints qStepNumber(iq) = 4 qLambda(iq) = 2 qMiu(iq) = 0 qinitValue(iq) = 10.0D0 qFinValue(iq) = 10.0D0 End Do numberConstraints_restart = 1 Allocate(qStepNumber_restart(1:numberConstraints_restart)) Allocate( qLambda_restart(1:numberConstraints_restart)) Allocate( qMiu_restart(1:numberConstraints_restart)) Allocate( qinitValue_restart(1:numberConstraints_restart)) Allocate( qFinValue_restart(1:numberConstraints_restart)) Do iq = 1, numberConstraints_restart qStepNumber_restart(iq) = 4 qLambda_restart(iq) = 2 qMiu_restart(iq) = 0 qinitValue_restart(iq) = 10.0D0 qFinValue_restart(iq) = 10.0D0 End Do neck_initValue = 1.0D0 neck_finValue = 1.0D0 neck_numberPoints = 1 neck_initValue_restart = 1.0D0 neck_finValue_restart = 1.0D0 neck_numberPoints_restart = 1 basis_frequency_initValue = 8.0D0 basis_frequency_stepValue = 0.1D0 basis_frequency_stepNumber = 1 basis_size_initValue = 8 basis_size_stepValue = 2 basis_size_stepNumber = 1 basis_states_initValue = 165 basis_states_stepValue = 2 basis_states_stepNumber = 1 basis_deform_initValue = 0.0D0 basis_deform_stepValue = 0.1D0 basis_deform_stepNumber = 1 End Subroutine mpi_initializeParallelData !-------------------------------------------------------------------! ! This subroutine sets variables that are used to construct the ! ! list of jobs. The subroutine is executed by all processors ! !-------------------------------------------------------------------! Subroutine mpi_setParallelData() Integer :: ii,jj,iq,iq_batch Character(Len=10) :: KEYWOR, string_10 !=================================================! ! ! ! RECONSTRUCTION OF THE INPUT OF TYPE INTEGER ! ! ! !=================================================! ii = 1 Do While ( ii .le. counterInt ) KEYWOR = vectorAll_names(ii) If(trim(KEYWOR).Eq.'CALCULMODE') Then modeMPI_basis = vectorInt_value(ii) ii = ii + 1 modeMPI_deformation = vectorInt_value(ii) End If If(trim(KEYWOR).Eq.'CALCULNECK') Then modeMPI_neck = vectorInt_value(ii) End If If(trim(KEYWOR).Eq.'CONSTR_DEF'.OR.trim(KEYWOR).Eq.'CONSTR-DEF') Then turnConstraintOn = vectorInt_value(ii) End If If(trim(KEYWOR).Eq.'BATCH_MODE'.OR.trim(KEYWOR).Eq.'BATCH-MODE') Then number_batch = vectorInt_value(ii) ii = ii + 1 batch_mode = vectorInt_value(ii) If(batch_mode.Eq.0) number_batch=1 End If If(trim(KEYWOR).Eq.'BATCH_SPEC'.OR.trim(KEYWOR).Eq.'BATCH-SPEC') Then number_constraint_batch = vectorInt_value(ii) If(Allocated(lambda_batch)) Then Deallocate(lambda_batch,miu_batch) End If Allocate(lambda_batch(1:number_constraint_batch)) Allocate( miu_batch(1:number_constraint_batch)) Do iq_batch = 1, number_constraint_batch ii = ii + 1 miu_batch(iq_batch) = vectorInt_value(ii) ii = ii + 1 lambda_batch(iq_batch) = Abs(vectorInt_value(ii)) End Do End If If(trim(KEYWOR).Eq.'ALL_NUCLEI'.OR.trim(KEYWOR).Eq.'ALL-NUCLEI') Then neutron_stepNumber = vectorInt_value(ii) ii = ii + 1 neutron_stepValue = vectorInt_value(ii) ii = ii + 1 neutron_initValue = vectorInt_value(ii) ii = ii + 1 proton_stepNumber = vectorInt_value(ii) ii = ii + 1 proton_stepValue = vectorInt_value(ii) ii = ii + 1 proton_initValue = vectorInt_value(ii) End If If(trim(KEYWOR).Eq.'REALNUCLEI') Then neutron_stepNumber = vectorInt_value(ii) ii = ii + 1 proton_stepNumber = vectorInt_value(ii) End If If(trim(KEYWOR).Eq.'MULTICONST') Then numberConstraints = vectorInt_value(ii) If(Allocated(qLambda)) Then Deallocate(qLambda,qMiu,qStepNumber) End If Allocate(qLambda(1:numberConstraints)) Allocate(qMiu(1:numberConstraints)) Allocate(qStepNumber(1:numberConstraints)) Do iq = 1, numberConstraints ii = ii + 1 qStepNumber(iq) = vectorInt_value(ii) ii = ii + 1 qMiu(iq) = vectorInt_value(ii) ii = ii + 1 qLambda(iq) = vectorInt_value(ii) End Do End If If(trim(KEYWOR).Eq.'MULTIRESTA') Then numberConstraints_restart = vectorInt_value(ii) If(Allocated(qLambda_restart)) Then Deallocate(qLambda_restart,qMiu_restart,qStepNumber_restart) End If Allocate( qLambda_restart(1:numberConstraints_restart)) Allocate( qMiu_restart(1:numberConstraints_restart)) Allocate(qStepNumber_restart(1:numberConstraints_restart)) Do iq = 1, numberConstraints_restart ii = ii + 1 qStepNumber_restart(iq) = vectorInt_value(ii) ii = ii + 1 qMiu_restart(iq) = vectorInt_value(ii) ii = ii + 1 qLambda_restart(iq) = vectorInt_value(ii) End Do End If If(trim(keywor) == 'NECK_CONST') Then neck_numberPoints = vectorInt_value(ii) End If If(trim(keywor) == 'NECK_RESTA') Then neck_numberPoints_restart = vectorInt_value(ii) End If If(trim(KEYWOR).Eq.'OPTIM_GRID'.OR.trim(KEYWOR).Eq.'OPTIM-GRID') Then optimizePath = vectorInt_value(ii) End If If(trim(KEYWOR).Eq.'ALL_FORCES'.OR.trim(KEYWOR).Eq.'ALL-FORCES') Then counterForce = vectorInt_value(ii) End If If(trim(KEYWOR).Eq.'TEMPRESTAR') Then restartTemperature = vectorInt_value(ii) End If If(trim(KEYWOR).Eq.'TEMPVALUES') Then counterTemperature = vectorInt_value(ii) End If If(trim(KEYWOR).Eq.'BASIS_FREQ'.OR.trim(KEYWOR).Eq.'BASIS-FREQ') Then basis_frequency_stepNumber = vectorInt_value(ii) End If If(trim(KEYWOR).Eq.'BASIS_NSHL'.OR.trim(KEYWOR).Eq.'BASIS-NSHL') Then basis_size_stepNumber = vectorInt_value(ii) ii = ii + 1 basis_size_stepValue = vectorInt_value(ii) ii = ii + 1 basis_size_initValue = vectorInt_value(ii) End If If(trim(KEYWOR).Eq.'BASIS_NSTA'.OR.trim(KEYWOR).Eq.'BASIS-NSTA') Then basis_states_stepNumber = vectorInt_value(ii) ii = ii + 1 basis_states_stepValue = vectorInt_value(ii) ii = ii + 1 basis_states_initValue = vectorInt_value(ii) End If If(trim(KEYWOR).Eq.'BASIS_DEFS'.OR.trim(KEYWOR).Eq.'BASIS-DEFS') Then basis_deform_stepNumber = vectorInt_value(ii) End If ! Go to the next item on the list ii = ii + 1 End Do Deallocate(vectorInt_value) !==============================================! ! ! ! RECONSTRUCTION OF THE INPUT OF TYPE REAL ! ! ! !==============================================! ii = 1 Do While ( ii .le. counterDbl ) KEYWOR = vectorAll_names(ii + counterInt) If(trim(KEYWOR).Eq.'REALNUCLEI') Then x_neutron_stepValue = vectorDbl_value(ii) ii = ii + 1 x_neutron_initValue = vectorDbl_value(ii) ii = ii + 1 x_proton_stepValue = vectorDbl_value(ii) ii = ii + 1 x_proton_initValue = vectorDbl_value(ii) End If If(trim(KEYWOR).Eq.'MULTICONST') Then If(Allocated(qinitValue)) Then Deallocate(qinitValue,qFinValue) End If Allocate(qinitValue(1:numberConstraints)) Allocate(qFinValue(1:numberConstraints)) Do iq = 1, numberConstraints qFinValue(iq) = vectorDbl_value(ii) ii = ii + 1 qinitValue(iq) = vectorDbl_value(ii) ii = ii + 1 End Do ii = ii - 1 End If If(trim(KEYWOR).Eq.'MULTIRESTA') Then If(Allocated(qinitValue_restart)) Then Deallocate(qinitValue_restart,qFinValue_restart) End If Allocate(qinitValue_restart(1:numberConstraints_restart)) Allocate( qFinValue_restart(1:numberConstraints_restart)) Do iq = 1, numberConstraints_restart qFinValue_restart(iq) = vectorDbl_value(ii) ii = ii + 1 qinitValue_restart(iq) = vectorDbl_value(ii) ii = ii + 1 End Do ii = ii - 1 End If If(trim(KEYWOR).Eq.'NECK_CONST') Then neck_finValue = vectorDbl_value(ii) ii = ii + 1 neck_initValue = vectorDbl_value(ii) End If If(trim(KEYWOR).Eq.'NECK_RESTA') Then neck_finValue_restart = vectorDbl_value(ii) ii = ii + 1 neck_initValue_restart = vectorDbl_value(ii) End If If(trim(KEYWOR).Eq.'TEMPVALUES') Then If(Allocated(temperatureVector)) Then Deallocate(temperatureVector) End If Allocate(temperatureVector(1:counterTemperature)) Do jj=1, counterTemperature temperatureVector(jj) = vectorDbl_value(ii+counterTemperature-jj) End Do ii = ii + counterTemperature - 1 End If If(trim(KEYWOR).Eq.'BASIS_FREQ'.OR.trim(KEYWOR).Eq.'BASIS-FREQ') Then basis_frequency_stepValue = vectorDbl_value(ii) ii = ii + 1 basis_frequency_initValue = vectorDbl_value(ii) End If If(trim(KEYWOR).Eq.'BASIS_DEFS'.OR.trim(KEYWOR).Eq.'BASIS-DEFS') Then basis_deform_stepValue = vectorDbl_value(ii) ii = ii + 1 basis_deform_initValue = vectorDbl_value(ii) End If ! Go to the next item on the list ii = ii + 1 End Do Deallocate(vectorDbl_value) !================================================! ! ! ! RECONSTRUCTION OF THE INPUT OF TYPE STRING*10 ! ! ! !================================================! ii = 1 Do While ( ii .le. counterStr ) KEYWOR = vectorAll_names(ii + counterInt + counterDbl) If(trim(KEYWOR).Eq.'ALL_FORCES'.OR.trim(KEYWOR).Eq.'ALL-FORCES') Then If(Allocated(forceVector)) Then Deallocate(forceVector) End If Allocate(forceVector(1:counterForce)) Do jj=1, counterForce string_10 = vectorStr_value(jj) forceVector(jj) = string_10(1:4) End Do End If ! Go to the next item on the list ii = ii + counterForce End Do Deallocate(vectorStr_value) Deallocate(vectorAll_names) End Subroutine mpi_setParallelData !-------------------------------------------------------------------! ! This subroutine prints variables that are common to all processes.! ! The subroutine is executed by processor 0 only ! !-------------------------------------------------------------------! Subroutine mpi_printSequentialData(IOPTGS,NXHERM,NYHERM,NZHERM, & IPAIRI,IPAHFB,IPABCS,IMFHFB, & LIPKIN,LIPKIP, & FE2FIN,IF2FIN,FE2FIP,IF2FIP, & IROTAT,IREVIE, & ICONTI,IPCONT,IGCONT,IGPCON, & ILCONT,IFCONT,IACONT, & NUCHAO,NUPING,EPSPNG,EPSITE, & ISIMPY,ISIGNY,IPARTY, & ISIMTX,ISIMTY,ISIMTZ, & NOITER,NULAST, & NOSCIL,NLIMIT,ENECUT, & ICOTYP,ICOUDI,ICOUEX, & I_GOGA,IGOGPA, & I1LINE, & SLOWEV,SLOWOD,SLOWPA,SLOWLI, & ECUTOF, & NILXYZ, & IBROYD,N_ITER,BROTRI,ALPHAM, & PRHO_T,PRHOST,POWERT, & PRHO_N,PRHOSN,PRHODN,POWERN, & PRHO_P,PRHOSP,PRHODP,POWERP, & OMEGAY, & IFSHEL,NPOLYN,GSTRUN,GSTRUP,HOMFAC, & LAMBAS,MIUBAS,BETBAS,NUMDEF, & LAM_WS,MIU_WS,NUM_WS,ITE_WS,STP_WS, & LMXINP,IFPUSH,IREAWS, & ISTAND,KETA_J,KETA_W,KETACM,KETA_M, & IF_RPA,IF_EDF,IF_THO,CBETHO, & IBASIS,IFRAGM, & NMUCON,NMUCOU,NMUPRI, & NEXBET,IPRIBE,IPRIBL, & FCHOM0, & IWRIFI,IFTEMP,IFNECK, & INNUMB,IZNUMB,R0PARM, & DELTAE,XLOCMX,V2_MIN,NTHETA,MIN_QP, & IDEALL,IDELOC,IDECON, & batch_init) Integer, intent(inout) :: ICOTYP, ICOUDI, ICOUEX, I_GOGA, IGOGPA, NOITER, IBROYD, & N_ITER, NULAST, NUPING, NUCHAO, ISIMPY, ISIGNY, IPARTY, & IROTAT, ISIMTX, ISIMTY, ISIMTZ, IPAIRI, IPAHFB, IPABCS, & IMFHFB, LIPKIN, LIPKIP, IOPTGS, NXHERM, NYHERM, NZHERM, & NOSCIL, NLIMIT, I1LINE, NILXYZ, IREVIE, ICONTI, IPCONT, & IGCONT, IGPCON, ILCONT, IFCONT, IACONT, IFSHEL, NPOLYN, & NUMDEF, NUM_WS, ITE_WS, IREAWS, ISTAND, KETA_J, KETA_W, & KETACM, KETA_M, IF_RPA, IF_EDF, IF_THO, IBASIS, IFRAGM, & IFPUSH, NMUCON, NMUCOU, NMUPRI, NEXBET, IPRIBE, IPRIBL, & IWRIFI, IFTEMP, IFNECK, LMXINP, INNUMB, IZNUMB, NTHETA, & MIN_QP, IDEALL, IDELOC, IDECON, batch_init,IF2FIN,IF2FIP Real, intent(inout) :: BROTRI, SLOWEV, SLOWOD, SLOWPA, SLOWLI, EPSITE, & EPSPNG, ECUTOF, ENECUT, PRHO_T, PRHOST, POWERT, & OMEGAY, GSTRUN, GSTRUP, HOMFAC, STP_WS, ALPHAM, & CBETHO, FCHOM0, PRHO_N, PRHOSN, PRHODN, POWERN, & PRHO_P, PRHOSP, PRHODP, POWERP, R0PARM, DELTAE, & XLOCMX, V2_MIN, FE2FIN, FE2FIP Integer :: ii,nfisum,LAMACT,MIUACT Integer :: mpi_rank,mpi_err Integer, Dimension(:), pointer :: LAMBAS,MIUBAS,LAM_WS,MIU_WS Real, Dimension(:), pointer :: BETBAS COMMON /OUTPUT/ nfisum Call mpi_comm_rank(MPI_COMM_WORLD,mpi_rank,mpi_err) ! Printing out the characteristics of the parallel run If(mpi_rank.Eq.0 .And. batch_init.Eq.1) Then Write(nfisum,'(/,14X,"/",42("-"),"/")') Write(nfisum,'(14X,"/ CHARACTERISTICS OF THE RUN: SERIAL PART /")') Write(nfisum,'(14X,"/",42("-"),"/",/)') Write(nfisum,'("/",77("-"),"/",/,"/",77X,"/")') Write(nfisum,'("/ Iterations",64X,"/")') Write(nfisum,'("/ ----------",64X,"/")') Write(nfisum,'("/ IBROYD = ",i1," N_ITER = ",i3," BROTRI = ",f9.3," ALPHAM = ",f7.4,13x, & & "/",/,"/",77X,"/")') IBROYD,N_ITER,BROTRI,ALPHAM Write(nfisum,'("/ NOITER = ",i5," NULAST = ",i5," NUCHAO = ",i5," NUPING = ",i5,13X,"/")') & NOITER,NULAST,NUCHAO,NUPING Write(nfisum,'("/ SLOWEV = ",f7.4," SLOWOD = ",f7.4," SLOWPA = ",f7.4," SLOWLI = ",f7.4,5x,"/")') & SLOWEV,SLOWOD,SLOWPA,SLOWLI Write(nfisum,'("/ EPSPNG = ",e20.12," EPSITE = ",e20.12,13X,"/")') & EPSPNG,EPSITE Write(nfisum,'("/",77X,"/")') Write(nfisum,'("/ Characteristics of the HO basis",43X,"/")') Write(nfisum,'("/ -------------------------------",43X,"/")') Write(nfisum,'("/ IBASIS = ",i3," IOPTGS = ",i1,49X,"/")') IBASIS,IOPTGS Write(nfisum,'("/ NOSCIL = ",i3," NLIMIT = ",i5," ENECUT = ",f9.3,26X,"/")') & NOSCIL,NLIMIT,ENECUT Write(nfisum,'("/ NXHERM = ",i3," NYHERM = ",i3," NZHERM = ",i3,34X,"/")') & NXHERM,NYHERM,NZHERM Write(nfisum,'("/ INNUMB = ",i3," IZNUMB = ",i5," R0PARM = ",f9.3,26X,"/")') & INNUMB,IZNUMB,R0PARM Write(nfisum,'("/",77X,"/")') If(Associated(LAMBAS)) Then Write(nfisum,'("/ Deformations of the nuclear surface",37X,"/")') Do ii=1,NUMDEF LAMACT=Abs(LAMBAS(ii)); MIUACT=MIUBAS(ii) Write(nfisum,'("/ l = ",i2," m = ",i2," beta_lm = ",f7.4,41X,"/")') & LAMACT,MIUACT,BETBAS(ii) End Do Write(nfisum,'("/",77X,"/")') End If Write(nfisum,'("/ Type of calculation",55X,"/")') Write(nfisum,'("/ -------------------",55X,"/")') Write(nfisum,'("/ IPAIRI = ",i1," IPAHFB = ",i1," IPABCS =",i2, & & " IMFHFB = ",i1," LIPKIN = ",i1," LIPKIP = ",i1,7X,"/")') & IPAIRI,IPAHFB,IPABCS,IMFHFB,LIPKIN,LIPKIP Write(nfisum,'("/ FE2FIN = ",f9.3," IF2FIN = ",i1," FE2FIP = ",f9.3," IF2FIN =",i1,7X,"/")') & FE2FIN,IF2FIN,FE2FIP,IF2FIP Write(nfisum,'("/ ICOUDI = ",i1," ICOUEX = ",i1," I_GOGA = ",i1, & & " IGOGPA = ",i1," IFTEMP = ",i1," IFSHEL = ",i1,7X,"/")') & ICOUDI,ICOUEX,I_GOGA,IGOGPA,IFTEMP,IFSHEL Write(nfisum,'("/ IFRAGM = ",i1," ICOUEX = ",i1,51X,"/")') & IFRAGM,IF_THO Write(nfisum,'("/ IF_EDF = ",i1," ISTAND = ",i1," KETA_J = ",i1, & & " KETA_W = ",i1," KETACM = ",i1," KETA_M = ",i1,7X,"/")') & IF_EDF,ISTAND,KETA_J,KETA_W,KETACM,KETA_M Write(nfisum,'("/",77X,"/")') Write(nfisum,'("/ Pairing options",59X,"/")') Write(nfisum,'("/ ---------------",59X,"/")') Write(nfisum,'("/ ECUTOF = ",f9.3,54X,"/")') ECUTOF Write(nfisum,'("/ PRHO_T = ",f9.3," PRHOST = ",f9.3," POWERT = ",f9.3,16X,"/")') & PRHO_T,PRHOST,POWERT Write(nfisum,'("/ PRHO_N = ",f9.3," PRHOSN = ",f9.3," POWERN = ",f9.3,16X,"/")') & PRHO_N,PRHOSN,POWERN Write(nfisum,'("/ PRHO_P = ",f9.3," PRHOSP = ",f9.3," POWERP = ",f9.3,16X,"/")') & PRHO_P,PRHOSP,POWERP Write(nfisum,'("/",77X,"/")') Write(nfisum,'("/ Constraints",63X,"/")') Write(nfisum,'("/ -----------",63X,"/")') Write(nfisum,'("/ IF_RPA = ",i2," IFPUSH = ",i2," IFNECK = ",i2,37X,"/")') & IF_RPA,IFPUSH,IFNECK Write(nfisum,'("/ NMUCON = ",i2," NMUCOU = ",i2," NMUPRI = ",i2,37X,"/")') & NMUCON,NMUCOU,NMUPRI Write(nfisum,'("/ OMEGAY = ",f7.4,56X,"/")') & OMEGAY Write(nfisum,'("/",77X,"/")') Write(nfisum,'("/ Restart options",59X,"/")') Write(nfisum,'("/ ---------------",59X,"/")') Write(nfisum,'("/ ICONTI = ",i1," IPCONT = ",i1," IGCONT = ",i1, & & " IGPCON = ",i1," ILCONT = ",i1," IFCONT = ",i1,7X,"/")') & ICONTI,IPCONT,IGCONT,IGPCON,ILCONT,IFCONT Write(nfisum,'("/ IREAWS = ",i1,62X,"/")') IREAWS If(IREAWS.Eq.1) Then Write(nfisum,'("/ Minimization over deformations",42X,"/")') Write(nfisum,'("/ NUMINP = ",i3," ITEINP = ",i3," LMXINP = ",i2," STPINP = ",f7.4,15X,"/")') & NUM_WS,ITE_WS,LMXINP,STP_WS Do ii=1,NUM_WS LAMACT=Abs(LAM_WS(ii)); MIUACT=MIU_WS(ii) Write(nfisum,'("/ l = ",i2," m = ",i2,59X,"/")') LAMACT,MIUACT End Do End If Write(nfisum,'("/",77X,"/")') Write(nfisum,'("/ Symmetry flags",60X,"/")') Write(nfisum,'("/ --------------",60X,"/")') Write(nfisum,'("/ IROTAT = ",i2," ISIMPY = ",i2," ISIGNY = ",i2," IPARTY = ",i2,25X,"/")') & IROTAT,ISIMPY,ISIGNY,IPARTY Write(nfisum,'("/ ISIMTX = ",i2," ISIMTY = ",i2," ISIMTZ = ",i2,25X,"/")') & ISIMTX,ISIMTY,ISIMTZ Write(nfisum,'("/",77X,"/")') Write(nfisum,'("/ Various",67X,"/")') Write(nfisum,'("/ -------",67X,"/")') Write(nfisum,'("/ I1LINE = ",i2," NILXYZ = ",i2," IREVIE = ",i2," IWRIFI = ",i2,25X,"/")') & I1LINE,NILXYZ,IREVIE,IWRIFI Write(nfisum,'("/ NEXBET = ",i2," IPRIBE = ",i2," IPRIBL = ",i2,37X,"/")') & NEXBET,IPRIBE,IPRIBL If(IFSHEL.Gt.0) then Write(nfisum,'("/",77X,"/")') Write(nfisum,'("/ NPOLYN = ",i2," HOMFAC = ",f7.4," GSTRUN = ",f7.4," GSTRUP = ",f7.4,10X,"/")') & NPOLYN,HOMFAC,GSTRUN,GSTRUP End if Write(nfisum,'("/ MIN_QP = ",i2," DELTAE = ",f7.3," XLOCMX = ",f7.5," V2_MIN = ",f7.5,10X,"/")') & MIN_QP,DELTAE,XLOCMX,V2_MIN Write(nfisum,'("/ IDEALL = ",i2," IDELOC = ",i2," IDECON = ",i2,37X,"/")') & IDEALL,IDELOC,IDECON Write(nfisum,'("/",77X,"/",/,"/",77("-"),"/")') End If Call mpi_barrier(MPI_COMM_WORLD, mpi_err) End Subroutine mpi_printSequentialData !-------------------------------------------------------------------! ! This subroutine prints variables that are used to construct the ! ! list of jobs. The subroutine is executed by processor 0 only ! !-------------------------------------------------------------------! Subroutine mpi_printParallelData() Integer :: ii,nfisum,mpi_rank,mpi_err COMMON /OUTPUT/ nfisum Call mpi_comm_rank(MPI_COMM_WORLD,mpi_rank,mpi_err) ! Printing out the characteristics of the parallel run If(mpi_rank.Eq.0) Then Write(nfisum,'(/,14X,"/",44("-"),"/")') Write(nfisum,'(14X,"/ CHARACTERISTICS OF THE RUN: PARALLEL PART /")') Write(nfisum,'(14X,"/",44("-"),"/",/)') Write(nfisum,'("/",77("-"),"/",/,"/",77X,"/")') Write(nfisum,'("/ Calculation mode",58X,"/")') Write(nfisum,'("/ ----------------",58X,"/")') Write(nfisum,'("/ Deformation .: ",i2,55X,"/")') modeMPI_deformation Write(nfisum,'("/ Neck ........: ",i2,55X,"/")') modeMPI_neck Write(nfisum,'("/ Batch .......: ",i2,55X,"/")') batch_mode Write(nfisum,'("/ Batch number : ",i2,55X,"/")') number_batch Write(nfisum,'("/ Basis .......: ",i2,55X,"/")') modeMPI_basis Write(nfisum,'("/",77X,"/")') Write(nfisum,'("/ Nucleus grid",62X,"/")') Write(nfisum,'("/ ------------",62X,"/")') Write(nfisum,'("/ Particle Initial Step Number of steps",28X,"/")') If(neutron_initValue.Gt.0.0D0) Then Write(nfisum,'("/ neutron ",F8.4," ",F7.4," ",i4,34X,"/")') & x_neutron_initValue, x_neutron_stepValue, neutron_stepNumber Else Write(nfisum,'("/ neutron ",i4," ",i2," ",i4,34X,"/")') & neutron_initValue, neutron_stepValue, neutron_stepNumber End If If(proton_initValue.Gt.0.0D0) Then Write(nfisum,'("/ proton ",F8.4," ",F7.4," ",i4,34X,"/")') & x_proton_initValue, x_proton_stepValue, proton_stepNumber Else Write(nfisum,'("/ proton ",i4," ",i2," ",i4,34X,"/")') & proton_initValue, proton_stepValue, proton_stepNumber End If Write(nfisum,'("/",77X,"/")') ! Displaying characteristics of the default deformation grid (regular) If(modeMPI_deformation.Eq.1.Or.modeMPI_basis.Eq.1) Then Write(nfisum,'("/ Deformation grid",58X,"/")') Write(nfisum,'("/ ----------------",58X,"/")') Write(nfisum,'("/ Constrained calculations: ",i4,42X,"/")') turnConstraintOn Write(nfisum,'("/ Number of constraints ..: ",i4,42X,"/")') numberConstraints Write(nfisum,'("/ Constraint # Lambda Miu Initial Final Number of steps",7X,"/")') Do ii=1,numberConstraints Write(nfisum,'("/ ",6X,i2,9X,i2,6X,i2,2X,f8.2,2X,f8.2,6X,i5,14X,"/")') & ii, qLambda(ii), qMiu(ii), qinitValue(ii), qFinValue(ii), qStepNumber(ii) End Do Write(nfisum,'("/",77X,"/")') End If ! Displaying characteristics of the initial grid (regular) and restart deformation grid (regular) If(modeMPI_deformation.Eq.2) Then Write(nfisum,'("/ Original deformation grid",49X,"/")') Write(nfisum,'("/ -------------------------",49X,"/")') Write(nfisum,'("/ Constrained calculations: ",i4,42X,"/")') turnConstraintOn Write(nfisum,'("/ Number of constraints ..: ",i4,42X,"/")') numberConstraints Write(nfisum,'("/ Constraint # Lambda Miu Initial Final Number of steps",7X,"/")') Do ii=1,numberConstraints Write(nfisum,'("/ ",6X,i2,9X,i2,6X,i2,2X,f8.2,2X,f8.2,6X,i5,14X,"/")') & ii, qLambda(ii), qMiu(ii), qinitValue(ii), qFinValue(ii), qStepNumber(ii) End Do Write(nfisum,'("/ New deformation grid",54X,"/")') Write(nfisum,'("/ --------------------",54X,"/")') Write(nfisum,'("/ Optimize path ..........: ",i4,42X,"/")') optimizePath Write(nfisum,'("/ Number of constraints ..: ",i4,42X,"/")') numberConstraints_restart Write(nfisum,'("/ Constraint # Lambda Miu Initial Final Number of steps",7X,"/")') Do ii=1,numberConstraints_restart Write(nfisum,'("/ ",6X,i2,9X,i2,6X,i2,2X,f8.2,2X,f8.2,6X,i5,14X,"/")') & ii, qLambda_restart(ii), qMiu_restart(ii), & qinitValue_restart(ii), qFinValue_restart(ii), qStepNumber_restart(ii) End Do Write(nfisum,'("/",77X,"/")') End If ! Displaying characteristics of the initial grid (only what constraints were used) and of the ! restart grid (regular) If(modeMPI_deformation.Eq.3) Then Write(nfisum,'("/ Original deformation grid",49X,"/")') Write(nfisum,'("/ -------------------------",49X,"/")') Write(nfisum,'("/ Number of constraints ..: ",i4,42X,"/")') numberConstraints Write(nfisum,'("/ Constraint # Lambda Miu",44X,"/")') Do ii=1,numberConstraints Write(nfisum,'("/ ",6X,i2,9X,i2,6X,i2,45X,"/")') & ii, qLambda(ii), qMiu(ii) End Do Write(nfisum,'("/ New deformation grid",54X,"/")') Write(nfisum,'("/ --------------------",54X,"/")') Write(nfisum,'("/ Optimize path ..........: ",i4,42X,"/")') optimizePath Write(nfisum,'("/ Number of constraints ..: ",i4,42X,"/")') numberConstraints_restart Write(nfisum,'("/ Constraint # Lambda Miu Initial Final Number of steps",7X,"/")') Do ii=1,numberConstraints_restart Write(nfisum,'("/ ",6X,i2,9X,i2,6X,i2,2X,f8.2,2X,f8.2,6X,i5,14X,"/")') & ii, qLambda_restart(ii), qMiu_restart(ii), & qinitValue_restart(ii), qFinValue_restart(ii), qStepNumber_restart(ii) End Do write(nfisum,'("/",77X,"/")') End If ! Displaying characteristics of the origina deformation grid (grid). New calculation will not ! yield a regular grid but an irregular path defined by what is in the file hfodd_path.d If(modeMPI_deformation.Eq.4) Then Write(nfisum,'("/ Original deformation grid",49X,"/")') Write(nfisum,'("/ -------------------------",49X,"/")') Write(nfisum,'("/ Number of constraints ..: ",i4,42X,"/")') numberConstraints Write(nfisum,'("/ Constraint # Lambda Miu",44X,"/")') Do ii=1,numberConstraints Write(nfisum,'("/ ",6X,i2,9X,i2,6X,i2,45X,"/")') & ii, qLambda(ii), qMiu(ii) End Do write(nfisum,'("/",77X,"/")') End If ! Displaying characteristics of the default neck grid If(modeMPI_neck.Eq.1) Then Write(nfisum,'("/ Neck grid",65X,"/")') Write(nfisum,'("/ ---------",65X,"/")') Write(nfisum,'("/ ",29X,f8.4,2X,f8.4,6X,i5,14X,"/")') & neck_initValue, neck_finValue, neck_numberPoints Write(nfisum,'("/",77X,"/")') End If ! Displaying characteristics of the initial and restart neck grid If(modeMPI_neck.Eq.2) Then Write(nfisum,'("/ Original neck grid",56X,"/")') Write(nfisum,'("/ ------------------",56X,"/")') Write(nfisum,'("/ ",29X,f8.4,2X,f8.4,6X,i5,14X,"/")') & neck_initValue, neck_finValue, neck_numberPoints Write(nfisum,'("/ New neck grid",61X,"/")') Write(nfisum,'("/ -------------",61X,"/")') Write(nfisum,'("/ ",29X,f8.4,2X,f8.4,6X,i5,14X,"/")') & neck_initValue_restart, neck_finValue_restart, neck_numberPoints_restart Write(nfisum,'("/",77X,"/")') End If Write(nfisum,'("/ Forces",68X,"/")') Write(nfisum,'("/ ------",68X,"/")') Do ii=1,counterForce Write(nfisum,'("/ ",A4,68X,"/")') forceVector(ii) End Do write(nfisum,'("/",77X,"/")') If(restartTemperature.Eq.1) Then Write(nfisum,'("/ Temperature Grid (restart is based on the exact same grid)",14X,"/")') Else Write(nfisum,'("/ Temperature Grid (restart is based on T=0 calculations)",18X,"/")') End If Write(nfisum,'("/ ----------------",58X,"/")') Do ii=1,counterTemperature Write(nfisum,'("/ ",F6.3,63X,"/")') temperatureVector(ii) End do Write(nfisum,'("/",77X,"/")') ! Displaying characteristics of the basis grid If(modeMPI_basis.Eq.1) Then Write(nfisum,'("/ Basis grid",64X,"/")') Write(nfisum,'("/ ----------",64X,"/")') Write(nfisum,'("/ Type Initial Step Number of steps",24X,"/")') Write(nfisum,'("/ Shells ",i2," ",i2," ",i2,32X,"/")') & basis_size_initValue,basis_size_stepValue,basis_size_stepNumber Write(nfisum,'("/ States ",i4," ",i4," ",i2,32X,"/")') & basis_states_initValue,basis_states_stepValue,basis_states_stepNumber Write(nfisum,'("/ HO frequency ",f6.3," ",f6.3," ",i2,32X,"/")') & basis_frequency_initValue,basis_frequency_stepValue,basis_frequency_stepNumber Write(nfisum,'("/ HO deformation ",f6.3," ",f6.3," ",i2,32X,"/")') & basis_deform_initValue,basis_deform_stepValue,basis_deform_stepNumber Write(nfisum,'("/",77X,"/")') End If Write(nfisum,'("/",77("-"),"/")') End If Call mpi_barrier(MPI_COMM_WORLD, mpi_err) End Subroutine mpi_printParallelData End Module hfodd_mpiio