Skip to content

Commit

Permalink
Merge pull request SPARC-X#206 from StevenZhangCSFM/NPT
Browse files Browse the repository at this point in the history
  • Loading branch information
phanish-suryanarayana authored Oct 31, 2023
2 parents 6f93914 + 21f65d6 commit 1b4c538
Show file tree
Hide file tree
Showing 37 changed files with 14,698 additions and 6,925 deletions.
6 changes: 6 additions & 0 deletions ChangeLog
Original file line number Diff line number Diff line change
Expand Up @@ -3,6 +3,12 @@
-Name
-changes

--------------
Oct 31, 2023
Name: Boqin Zhang
Changes: (md.c)
1. Add NPT_SCALE_CONSTRAINTS

--------------
Oct 28, 2023
Name: Xin Jing
Expand Down
1 change: 1 addition & 0 deletions doc/.LaTeX/Introduction.tex
Original file line number Diff line number Diff line change
Expand Up @@ -346,6 +346,7 @@
\hyperlink{NPT_NP_QMASS}{\texttt{NPT\_NP\_QMASS}} $\vert$
\hyperlink{NPT_NP_BMASS}{\texttt{NPT\_NP\_BMASS}} $\vert$
\hyperlink{NPT_SCALE_VECS}{\texttt{NPT\_SCALE\_VECS}} $\vert$
\hyperlink{NPT_SCALE_CONSTRAINTS}{\texttt{NPT\_SCALE\_CONSTRAINTS}} $\vert$
\hyperlink{TARGET_PRESSURE}{\texttt{TARGET\_PRESSURE}} $\vert$
\hyperlink{RESTART_FLAG}{\texttt{RESTART\_FLAG}} $\vert$
\hyperlink{TWTIME}{\texttt{TWTIME}}
Expand Down
52 changes: 48 additions & 4 deletions doc/.LaTeX/MD.tex
Original file line number Diff line number Diff line change
Expand Up @@ -72,7 +72,7 @@
\end{columns}

\begin{block}{Description}
Type of QMD to be performed. Currently, NVE (microcanonical ensemble), NVT\_NH (canonical ensemble with Nose-Hoover thermostat), NVK\_G (isokinetic ensemble with Gaussian thermostat), NPT\_NH (isothermal-isobaric ensemble with Nose-Hoover thermostat) and NPT\_NP (isothermal-isobaric ensemble with Nose-Poincare thermostat) are supported
Type of QMD to be performed. Currently, NVE (microcanonical ensemble), NVT\_NH (canonical ensemble with Nose-Hoover thermostat), NVK\_G (isokinetic ensemble with Gaussian thermostat), NPT\_NH (isothermal-isobaric ensemble with Nose-Hoover thermostat) and NPT\_NP (isothermal-isobaric ensemble with Nose-Poincare thermostat) are supported.
\end{block}

\end{frame}
Expand Down Expand Up @@ -517,7 +517,7 @@
\begin{columns}
\column{0.4\linewidth}
\begin{block}{Type}
Permutation of 1, 2, 3
Int
\end{block}

\begin{block}{Default}
Expand All @@ -535,11 +535,55 @@
\end{columns}

\begin{block}{Description}
Specify which lattice vectors can be rescaled in NPT-NH
Specify which lattice vectors can be rescaled in NPT\_NH and NPT\_NP. The cell will only expand or shrink in the specified directions.

Rescaled vectors can be specified for orthogonal systems if NPT\_NP thermostat is used.
\end{block}

\begin{block}{Remark}
Only numbers 1, 2 and 3 can be accepted. For example, if `` 2 3" is the input, the cell will only expand or shrink in the directions of lattice vector 2 and lattice vector 3. Only NPT-NH can specify the rescaled vector
Only three numbers 1, 2 and 3 can be accepted. For example, if ``2 3'' is the input, the cell will only expand or shrink in the directions of lattice vector 2 and lattice vector 3.

If it is set in NPT\_NH, the expansion or shrinkage on designated lattice vector will try to keep the total pressure to oscillate near the target pressure.

If it is set in NPT\_NP, the expansion or shrinkage on designated lattice vector will only try to keep the normal stress at their direction to oscillate near the target pressure.
\end{block}

\end{frame}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%



%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\begin{frame}[allowframebreaks]{\texttt{NPT\_SCALE\_CONSTRAINTS}} \label{NPT_SCALE_CONSTRAINTS}
\vspace*{-12pt}
\begin{columns}
\column{0.4\linewidth}
\begin{block}{Type}
Double
\end{block}

\begin{block}{Default}
none
\end{block}

\column{0.4\linewidth}
\begin{block}{Unit}
No unit
\end{block}

\begin{block}{Example}
\texttt{NPT\_SCALE\_CONSTRAINTS}: 12
\end{block}
\end{columns}

\begin{block}{Description}
Set the scale constraint for lattice vectors in NPT\_NP. The length ratio between the designated lattice vector keeps constant in NPT\_NP thermostat. For example, if ``12'' is set, then the length ratio between 1st and 2nd lattice vectors will keep constant.
\end{block}

\begin{block}{Remark}
Applicable to orthogonal system using NPT\_NP \hyperlink{MD_METHOD}{\texttt{MD\_METHOD}} only.

There are 4 types of available constraints. ``12'' or ``21''; ``13'' or ``31''; ``23'' or ``32''; ``123'' or ``132'' or ``213'' or ``231'' or ``312'' or ``321''.
\end{block}

\end{frame}
Expand Down
Binary file modified doc/Manual.pdf
Binary file not shown.
8 changes: 7 additions & 1 deletion src/include/isddft.h
Original file line number Diff line number Diff line change
Expand Up @@ -680,7 +680,11 @@ typedef struct _SPARC_OBJ{
double fs2atu; // conversion factor for femto second -> atomic unit of time (Jiffy)
// NPT
int NPTscaleVecs[3]; // which lattice vector can be rescaled?
int NPTisotropicFlag; // whether it is an isotropic cell expansion
int NPTconstraintFlag; // confinement on side length of cell. none: no length confinement (default); 1: a:b keeps unchanged; 2: a:c keeps unchanged;
// 3: a:c keeps unchanged; 4: a:b:c keeps unchanged, isotropic expansion. It is only available for NPT_NP.
int NPTisotropicFlag; // whether it is an isotropic cell expansion; a:b:c keeps similar during NPT.
// For NPT_NH, if all 3 lattive vectors are scalable, it will be an isotropic expansion;
// For NPT_NP, if all 3 lattive vectors are scalable, AND NPTconstraintFlag is 4, it will be an isotropic expansion.
double prtarget; // Target pressure of barostatic system, used in both NPT_NH and NPT_NP
double scale; // length ratio of the size of cell in NPT, used in both NPT_NH and NPT_NP
double volumeCell; // volume of the cell, used in both NPT_NH and NPT_NP
Expand Down Expand Up @@ -1167,6 +1171,8 @@ typedef struct _SPARC_INPUT_OBJ{

int NPT_NHnnos; // number of thermostat variables in NPT_NH
int NPTscaleVecs[3]; // which lattice vector can be rescaled?
int NPTconstraintFlag; // confinement on side length of cell. none: no length confinement (default); 1: a:b keeps unchanged; 2: a:c keeps unchanged;
// 3: a:c keeps unchanged; 4: a:b:c keeps unchanged, isotropic expansion. It is only available for NPT_NP.
double NPT_NHqmass[L_QMASS];// qmass used in NPT_NH
double NPT_NHbmass; // Bmass used in NPT_NH
double prtarget; // Target pressure of barostatic system, UNIT on input file is GPa
Expand Down
92 changes: 78 additions & 14 deletions src/initialization.c
Original file line number Diff line number Diff line change
Expand Up @@ -44,7 +44,7 @@
#define min(x,y) ((x)<(y)?(x):(y))
#define max(x,y) ((x)>(y)?(x):(y))

#define N_MEMBR 174
#define N_MEMBR 175


/**
Expand Down Expand Up @@ -691,6 +691,7 @@ void set_defaults(SPARC_INPUT_OBJ *pSPARC_Input, SPARC_OBJ *pSPARC) {
pSPARC_Input->NPTscaleVecs[0] = 1;
pSPARC_Input->NPTscaleVecs[1] = 1;
pSPARC_Input->NPTscaleVecs[2] = 1; // default lattice vectors to be rescaled in NPT
pSPARC_Input->NPTconstraintFlag = 0; // confinement on side length of cell. none: no length confinement (default)
pSPARC_Input->NPT_NHnnos = 0; // default amount of thermo variable for NPT_NH. If MDMeth is this but nnos is 0, program will stop
for (int subscript_NPTNH_qmass = 0; subscript_NPTNH_qmass < L_QMASS; subscript_NPTNH_qmass++){
pSPARC_Input->NPT_NHqmass[subscript_NPTNH_qmass] = 0.0;
Expand Down Expand Up @@ -1211,6 +1212,7 @@ void SPARC_copy_input(SPARC_OBJ *pSPARC, SPARC_INPUT_OBJ *pSPARC_Input) {
pSPARC->NPTscaleVecs[0] = pSPARC_Input->NPTscaleVecs[0];
pSPARC->NPTscaleVecs[1] = pSPARC_Input->NPTscaleVecs[1];
pSPARC->NPTscaleVecs[2] = pSPARC_Input->NPTscaleVecs[2];
pSPARC->NPTconstraintFlag = pSPARC_Input->NPTconstraintFlag;
pSPARC->NPT_NHnnos = pSPARC_Input->NPT_NHnnos;
pSPARC->ion_elec_eqT = pSPARC_Input->ion_elec_eqT;
pSPARC->ion_vel_dstr = pSPARC_Input->ion_vel_dstr;
Expand Down Expand Up @@ -1354,16 +1356,6 @@ void SPARC_copy_input(SPARC_OBJ *pSPARC, SPARC_INPUT_OBJ *pSPARC_Input) {
// find exchange correltaion decomposition
xc_decomposition(pSPARC);

// check MDMeth availability
if ((strcmpi(pSPARC->MDMeth,"NVT_NH") && strcmpi(pSPARC->MDMeth,"NVE")
&& strcmpi(pSPARC->MDMeth,"NVK_G") && strcmpi(pSPARC->MDMeth,"NPT_NH") && strcmpi(pSPARC->MDMeth,"NPT_NP")) != 0) {
if (!rank){
printf("\nCannot recognize MDMeth = \"%s\"\n",pSPARC->MDMeth);
printf("MDMeth (MD Method) must be one of the following:\n\tNVT_NH\t NVE\t NVK_G\t NPT_NH\t NPT_NP\n");
}
exit(EXIT_FAILURE);
}

/* process the data read from input files */
Ntypes = pSPARC->Ntypes;

Expand Down Expand Up @@ -2570,6 +2562,65 @@ void SPARC_copy_input(SPARC_OBJ *pSPARC, SPARC_INPUT_OBJ *pSPARC_Input) {
pSPARC->PrintPsiFlag[3] = 0; pSPARC->PrintPsiFlag[4] = pSPARC->Nkpts-1; // k-point start/end index
pSPARC->PrintPsiFlag[5] = 0; pSPARC->PrintPsiFlag[6] = pSPARC->Nstates-1; // band start/end index
}

// check MDMeth availability
if ((strcmpi(pSPARC->MDMeth,"NVT_NH") && strcmpi(pSPARC->MDMeth,"NVE")
&& strcmpi(pSPARC->MDMeth,"NVK_G") && strcmpi(pSPARC->MDMeth,"NPT_NH") && strcmpi(pSPARC->MDMeth,"NPT_NP")) != 0) {
if (!rank){
printf("\nCannot recognize MDMeth = \"%s\"\n",pSPARC->MDMeth);
printf("MDMeth (MD Method) must be one of the following:\n\tNVT_NH\t NVE\t NVK_G\t NPT_NH\t NPT_NP\n");
}
exit(EXIT_FAILURE);
}
if (strcmpi(pSPARC->MDMeth,"NPT_NP") == 0) {
if (pSPARC->cell_typ > 10 && pSPARC->cell_typ < 20) { // check conflict for non-orthogonal cell systems
if (! (pSPARC->NPTscaleVecs[0] * pSPARC->NPTscaleVecs[1] * pSPARC->NPTscaleVecs[2])) {
if (!rank) {
printf("\nCurrently NPT_NP only support isotropic expansion for non-orthogonal cells. Please set NPT_SCALE_VECS: 1 2 3 \n");
printf("then set NPT_SCALE_CONSTRAINTS: 123 \n");
}
exit(EXIT_FAILURE);
}
if (pSPARC->NPTconstraintFlag != 4) {
if (!rank) {
printf("\nCurrently NPT_NP only support isotropic expansion for non-orthogonal cells. Please add or change NPT_SCALE_CONSTRAINTS: 123");
}
exit(EXIT_FAILURE);
}
}
if (pSPARC->NPTconstraintFlag == 1) { // check conflict between NPT_SCALE_CONFINEMENTS and NPT_SCALE_VECS
if (! (pSPARC->NPTscaleVecs[0] * pSPARC->NPTscaleVecs[1])) { // a or b cannot be rescaled
if (!rank) {
printf("\nNPT_SCALE_CONSTRAINTS 12 has conflict with NPT_SCALE_VECS!\n");
}
exit(EXIT_FAILURE);
}
}
if (pSPARC->NPTconstraintFlag == 2) {
if (! (pSPARC->NPTscaleVecs[0] * pSPARC->NPTscaleVecs[2])) { // a or c cannot be rescaled
if (!rank) {
printf("\nNPT_SCALE_CONSTRAINTS 13 has conflict with NPT_SCALE_VECS!\n");
}
exit(EXIT_FAILURE);
}
}
if (pSPARC->NPTconstraintFlag == 3) {
if (! (pSPARC->NPTscaleVecs[1] * pSPARC->NPTscaleVecs[2])) { // b or c cannot be rescaled
if (!rank) {
printf("\nNPT_SCALE_CONSTRAINTS 23 has conflict with NPT_SCALE_VECS!\n");
}
exit(EXIT_FAILURE);
}
}
if (pSPARC->NPTconstraintFlag == 4) {
if (! (pSPARC->NPTscaleVecs[0] * pSPARC->NPTscaleVecs[1] * pSPARC->NPTscaleVecs[2])) { // a or b or c cannot be rescaled
if (!rank) {
printf("\nNPT_SCALE_CONSTRAINTS 123 has conflict with NPT_SCALE_VECS!\n");
}
exit(EXIT_FAILURE);
}
}
}
}


Expand Down Expand Up @@ -3338,7 +3389,7 @@ void write_output_init(SPARC_OBJ *pSPARC) {
}

fprintf(output_fp,"***************************************************************************\n");
fprintf(output_fp,"* SPARC (version Oct 28, 2023) *\n");
fprintf(output_fp,"* SPARC (version Oct 31, 2023) *\n");
fprintf(output_fp,"* Copyright (c) 2020 Material Physics & Mechanics Group, Georgia Tech *\n");
fprintf(output_fp,"* Distributed under GNU General Public License 3 (GPL) *\n");
fprintf(output_fp,"* Start time: %s *\n",c_time_str);
Expand Down Expand Up @@ -3482,6 +3533,17 @@ void write_output_init(SPARC_OBJ *pSPARC) {
}
if(strcmpi(pSPARC->MDMeth,"NPT_NP") == 0) {
//fprintf(output_fp,"AMOUNT_THERMO_VARIABLE: %d\n",pSPARC->NPT_NHnnos);
fprintf(output_fp,"NPT_SCALE_VECS:");
if (pSPARC->NPTscaleVecs[0] == 1) fprintf(output_fp," 1");
if (pSPARC->NPTscaleVecs[1] == 1) fprintf(output_fp," 2");
if (pSPARC->NPTscaleVecs[2] == 1) fprintf(output_fp," 3");
fprintf(output_fp,"\n");
fprintf(output_fp,"NPT_SCALE_CONSTRAINTS:");
if (pSPARC->NPTconstraintFlag == 0) fprintf(output_fp," none\n");
else if (pSPARC->NPTconstraintFlag == 1) fprintf(output_fp," 12\n");
else if (pSPARC->NPTconstraintFlag == 2) fprintf(output_fp," 13\n");
else if (pSPARC->NPTconstraintFlag == 3) fprintf(output_fp," 23\n");
else if (pSPARC->NPTconstraintFlag == 4) fprintf(output_fp," 123\n");
fprintf(output_fp,"NPT_NP_QMASS: %.15g\n",pSPARC->NPT_NP_qmass);
fprintf(output_fp,"NPT_NP_BMASS: %.15g\n",pSPARC->NPT_NP_bmass);
fprintf(output_fp,"TARGET_PRESSURE: %.15g GPa\n",pSPARC->prtarget);
Expand Down Expand Up @@ -3849,7 +3911,7 @@ void SPARC_Input_MPI_create(MPI_Datatype *pSPARC_INPUT_MPI) {
MPI_INT, MPI_INT, MPI_INT, MPI_INT, MPI_INT,
MPI_INT, MPI_INT, MPI_INT, MPI_INT, MPI_INT,
MPI_INT, MPI_INT, MPI_INT, MPI_INT, MPI_INT,
MPI_INT, MPI_INT, MPI_INT,
MPI_INT, MPI_INT, MPI_INT, MPI_INT,
MPI_DOUBLE, MPI_DOUBLE, MPI_DOUBLE, MPI_DOUBLE, MPI_DOUBLE,
MPI_DOUBLE,
MPI_DOUBLE, MPI_DOUBLE, MPI_DOUBLE, MPI_DOUBLE, MPI_DOUBLE,
Expand Down Expand Up @@ -3885,7 +3947,8 @@ void SPARC_Input_MPI_create(MPI_Datatype *pSPARC_INPUT_MPI) {
1, 1, 1, 1, 1,
1, 1, 1, 1, 1,
1, 1, 1, 1, 1,
1, 1, 1, 1, 1, 1, 1, 1, /* int */
1, 1, 1, 1, 1,
1, 1, 1, 1,/* int */
9, 3, L_QMASS, L_kpoint, L_kpoint,
L_kpoint, /* double array */
1, 1, 1, 1, 1,
Expand Down Expand Up @@ -3990,6 +4053,7 @@ void SPARC_Input_MPI_create(MPI_Datatype *pSPARC_INPUT_MPI) {
MPI_Get_address(&sparc_input_tmp.Poisson_solver, addr + i++);
MPI_Get_address(&sparc_input_tmp.d3Flag, addr + i++);
MPI_Get_address(&sparc_input_tmp.NPT_NHnnos, addr + i++);
MPI_Get_address(&sparc_input_tmp.NPTconstraintFlag, addr + i++);
MPI_Get_address(&sparc_input_tmp.MAXIT_FOCK, addr + i++);
MPI_Get_address(&sparc_input_tmp.EXXMeth_Flag, addr + i++);
MPI_Get_address(&sparc_input_tmp.ACEFlag, addr + i++);
Expand Down
Loading

0 comments on commit 1b4c538

Please sign in to comment.