Skip to content
Draft
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
11 changes: 6 additions & 5 deletions Makefile
Original file line number Diff line number Diff line change
@@ -1,14 +1,12 @@
GF90 = gfortran

XFLAG = -fopenmp -ffree-line-length-none -O3 -march=native -fimplicit-none -fwhole-file

LIBS = -llapack
LIBS = -llapack -lblas

BIN = Thomson

FFFLAGS = -Wall -g -msse4.2 -fcheck=all -Waliasing -Wampersand -Wconversion -Wsurprising -Wintrinsics-std -Wno-tabs -Wintrinsic-shadow -Wline-truncation -Wreal-q-constant

FFLAGS = -Wall -Wno-unused -Wno-unused-dummy-argument -O2
FFLAGS = -Wall -Wno-unused -Wno-unused-dummy-argument -O3 -fallow-argument-mismatch -g -mcpu=native -funroll-loops

help:
@echo " -----------------------------------------\n to make the code you should do : \n ----------------------------------------- \n make Thomson | to compiled with gfortran \n -----------------------------------------";\
Expand All @@ -22,4 +20,7 @@ clean:
$(RM) $(BIN) animation.mp4 data.out figure.png test fort.4 fort.3 energy.dat data_frame.dat

Thomson_d:
$(GF90) src/*.f90 -o $(BIN) $(FFFLAGS) $(LIBS)
$(GF90) src/*.f90 -o $(BIN) $(FFFLAGS) $(LIBS)

Thomson_mpi:
mpif90 src/*.f90 -o $(BIN) $(FFLAGS) $(LIBS)
176 changes: 112 additions & 64 deletions src/Thomson.f90
Original file line number Diff line number Diff line change
Expand Up @@ -2,6 +2,8 @@ program Thomson

implicit none

include 'mpif.h'

! ----------- ! variable declaration ! ----------- !

! ---- i ---- !
Expand Down Expand Up @@ -32,16 +34,32 @@ program Thomson
character (len = 20 ) :: show , hess , distance , animation
character (len = 20 ) :: origin
!-----------------------------------------------------------------------------------!

!-----------------------------------------------------------------------------------!
! MPI !

integer :: ierr
integer :: nprocs ! total number of processes
integer :: ME ! CPU index
character(len=2) :: Ncore


call MPI_Init(ierr)

call MPI_Comm_size(MPI_COMM_WORLD,nprocs,ierr) ! Nombre of procs dispo

call MPI_Comm_rank(MPI_COMM_WORLD,ME,ierr) ! Id de chaque proc


!-----------------------------------------------------------------------------------!

if (ME == 0) then

call system('clear')
if(command_argument_count().eq.0) then
print *, ""
write(*,'(a)')
write(*,'(a)') "Sorry, you didn't add the input file"
write(*,'(a)') "----------------------------------------------------"
write(*,'(a)') "you can run the program using ./Thomson << the name of the input file >>"
print *, ""
write(*,'(a)')
stop
else
call get_command_argument(1,arg)
Expand All @@ -51,10 +69,26 @@ program Thomson

call read_f(arg,space,n_ele,tol,itermax,Lx,Ly,Lz,typ,multi,show,hess,distance,animation,origin,nlines)

! ----- The Title ----- !
end if

call title(n_ele,space,itermax,typ,Lx,Ly,Lz)

call MPI_Bcast(space ,1,mpi_INTEGER,0,MPI_COMM_WORLD,ierr)
call MPI_Bcast(n_ele ,1,mpi_INTEGER,0,MPI_COMM_WORLD,ierr)
call MPI_Bcast(itermax,1,mpi_INTEGER,0,MPI_COMM_WORLD,ierr)
call MPI_Bcast(nlines ,1,mpi_INTEGER,0,MPI_COMM_WORLD,ierr)

call MPI_Bcast(tol ,1,mpi_double_precision,0,MPI_COMM_WORLD,ierr)
call MPI_Bcast(Lx ,1,mpi_double_precision,0,MPI_COMM_WORLD,ierr)
call MPI_Bcast(Ly ,1,mpi_double_precision,0,MPI_COMM_WORLD,ierr)
call MPI_Bcast(Lz ,1,mpi_double_precision,0,MPI_COMM_WORLD,ierr)

call MPI_Bcast(arg ,LEN(arg),mpi_CHARACTER,0,MPI_COMM_WORLD,ierr)
call MPI_Bcast(typ ,LEN(arg),mpi_CHARACTER,0,MPI_COMM_WORLD,ierr)
call MPI_Bcast(multi ,LEN(arg),mpi_CHARACTER,0,MPI_COMM_WORLD,ierr)
call MPI_Bcast(show ,LEN(arg),mpi_CHARACTER,0,MPI_COMM_WORLD,ierr)
call MPI_Bcast(hess ,LEN(arg),mpi_CHARACTER,0,MPI_COMM_WORLD,ierr)
call MPI_Bcast(distance ,LEN(arg),mpi_CHARACTER,0,MPI_COMM_WORLD,ierr)
call MPI_Bcast(animation,LEN(arg),mpi_CHARACTER,0,MPI_COMM_WORLD,ierr)
call MPI_Bcast(origin ,LEN(arg),mpi_CHARACTER,0,MPI_COMM_WORLD,ierr)

! ---- ! memory allocation ! --- !

Expand All @@ -78,11 +112,19 @@ program Thomson

call first_geo(arg,n_ele,space,Lx,Ly,Lz,typ,multi,geo,nlines)

write(*,*) ""
write(*,'(a)') " (( The starting Geometry )) "
write(*,*) ""
write(Ncore,'(I2)') ME

call print_geo(n_ele,space,geo)
open(10,file='rand_'//Ncore//'.out')

! ----- The Title ----- !

call title(n_ele,space,itermax,typ,Lx,Ly,Lz,10)

write(10,*) ""
write(10,'(a)') " (( The starting Geometry )) "
write(10,*) ""

call print_geo(n_ele,space,geo,10)

call energy(n_ele,space,geo,Lx,Ly,Lz,E)

Expand All @@ -93,21 +135,21 @@ program Thomson
! --------- ! condition if NAN or fall ! -------- !

if (isnan(Norm2(dervative_b)) .or. Norm2(dervative_b) > 1e20) then
write(*,'(a)') "____________________________________________________________________________________________"
write(*,'(a)') ""
write (*,'(a)') "you have 2 or more electron at the same position or the energy is infinity"
write(10,'(a)') "____________________________________________________________________________________________"
write(10,'(a)') ""
write(10,'(a)') "you have 2 or more electron at the same position or the energy is infinity"
Stop
end if

! --------- ! ! ----------- !

write(*,'(a)') "____________________________________________________________________________________________"
write(*,'(a)') ""
write(*,'(a,f24.16)') "The energy before optimization = " , E
write(*,'(a)') ""
write(*,'(a,f24.16)') "and the gradiant normalization = " , norm2(dervative_b)
write(*,'(a)') "____________________________________________________________________________________________"
write(*,'(a)') ""
write(10,'(a)') "____________________________________________________________________________________________"
write(10,'(a)') ""
write(10,'(a,f24.16)') "The energy before optimization = " , E
write(10,'(a)') ""
write(10,'(a,f24.16)') "and the gradiant normalization = " , norm2(dervative_b)
write(10,'(a)') "____________________________________________________________________________________________"
write(10,'(a)') ""

! ----- ! condition for only one calculation ! ---- !

Expand All @@ -120,15 +162,17 @@ program Thomson

if (show /= "show" .and. norm2(dervative_b) > tol) then

write(*,'(a)') ' Loop | Energy | Gradient NORM '
write(*,'(a)') "_______________________________________________________________"
write(*,'(a)') ""
write(10,'(a)') ' Loop | Energy | Gradient NORM '
write(10,'(a)') "_______________________________________________________________"
write(10,'(a)') ""
end if

if (animation == "animation") then
open(4,file='data_frame.dat',status = 'replace')
open(8,file='energy.dat' ,status = 'replace')
open(4,file='data_frame'//Ncore//'.dat',status = 'replace')
open(8,file='energy'//Ncore//'.dat' ,status = 'replace')
end if


! %%%%%%%%%%%%%%%%%%%%%%%%%%% !
! ---- ! the main loop ! ---- !
! %%%%%%%%%%%%%%%%%%%%%%%%%%% !
Expand All @@ -137,19 +181,18 @@ program Thomson

iter = iter + 1


if (norm2(dervative_b) < tol) then

call energy(n_ele,space,geo,Lx,Ly,Lz,E)

call hessian(n_ele,space,geo,Lx,Ly,Lz,H)

write(*,'(a)') ""
write(*,'(a,f24.16)') "you are at the minimum, The energy = " , E
write(*,'(a)') ""
write(*,'(a,E20.10,a)') "and the gradiant normalization = " , norm2(dervative_b), " smaller than the tolerance"
write(*,'(a)') "____________________________________________________________________________________________"
write(*,'(a)') ""
write(10,'(a)') ""
write(10,'(a,f24.16)') "you are at the minimum, The energy = " , E
write(10,'(a)') ""
write(10,'(a,E20.10,a)') "and the gradiant normalization = " , norm2(dervative_b), " smaller than the tolerance"
write(10,'(a)') "____________________________________________________________________________________________"
write(10,'(a)') ""
if (animation == "animation") then
call anim(n_ele,space,geo,E,iter,origin,Lx,Ly,Lz)
end if
Expand All @@ -164,13 +207,13 @@ program Thomson

if (norm2(dervative_a) < tol) then

write(*,'(a)') ""
write(*,'(a)') "____________________________________________________________________________________________"
write(*,'(a)') ""
write(*,'(a,f24.16,a,i5,a)') "Geometry converged at ", E ," after " ,iter-1 , " Loops"
write(*,'(a)') ""
write(*,'(a,E16.10)') "The gradiant norm = ", norm2(dervative_a)
write(*,'(a)') ""
write(10,'(a)') ""
write(10,'(a)') "____________________________________________________________________________________________"
write(10,'(a)') ""
write(10,'(a,f24.16,a,i5,a)') "Geometry converged at ", E ," after " ,iter-1 , " Loops"
write(10,'(a)') ""
write(10,'(a,E16.10)') "The gradiant norm = ", norm2(dervative_a)
write(10,'(a)') ""
exit
end if

Expand Down Expand Up @@ -208,21 +251,21 @@ program Thomson

if (show == "show") then

write(*,'(a,i7)') "Loop = " , iter
write(*,'(a)') ""
write(*,'(a)') "The geometry after the optimization:"
write(10,'(a,i7)') "Loop = " , iter
write(10,'(a)') ""
write(10,'(a)') "The geometry after the optimization:"

call print_geo(n_ele,space,geo)
call print_geo(n_ele,space,geo,10)

write(*,'(a)') "____________________________________________________________________________________________"
write(*,'(a)') ""
write(*,'(a,f24.16,a,E8.2)') "The energy after optimization = ", E , " | The derivative = ", norm2(dervative_a)
write(*,'(a)') "____________________________________________________________________________________________"
write(*,'(a)') ""
write(10,'(a)') "____________________________________________________________________________________________"
write(10,'(a)') ""
write(10,'(a,f24.16,a,E8.2)') "The energy after optimization = ", E , " | The derivative = ", norm2(dervative_a)
write(10,'(a)') "____________________________________________________________________________________________"
write(10,'(a)') ""

else

write(*,'(i7,(1x,f24.16),(3x,f24.16))') iter , E , norm2(dervative_a)
write(10,'(i7,(1x,f24.16),(3x,f24.16))') iter , E , norm2(dervative_a)

end if

Expand All @@ -234,39 +277,40 @@ program Thomson

if (show /= "show") then

write(*,'(a)') ""
write(*,'(a)') "The final geometry after the optimization:"
call print_geo(n_ele,space,geo)
write(10,'(a)') ""
write(10,'(a)') "The final geometry after the optimization:"

call print_geo(n_ele,space,geo,10)

end if


if (hess == "hessian") then

write(*,'(a)') "____________________________________________________________________________________________"
write(*,'(a)') ""
write(*,'(a)') " The Hessian matrix after the converge "
write(*,'(a)') ""

do i=1,n_ele*space
write (*,'(100(3x,f16.8))') H(i,:)
end do
write(10,'(a)') "____________________________________________________________________________________________"
write(10,'(a)') ""
write(10,'(a)') " The Hessian matrix after the converge "
write(10,'(a)') ""

call matout(n_ele*space,n_ele*space,H,10)

end if

if (distance == "distance") then
call print_distance(n_ele,space,geo,Lx,Ly,Lz)
call print_distance(n_ele,space,geo,Lx,Ly,Lz,10)
end if

call diagonalize_matrix(n_ele*space,H,Eign)
call diagonalize_matrix(n_ele*space,H,Eign,10)

if (animation == "animation") then
close(4)
close(8)
end if

if (animation == "animation") then
call plot_anim(n_ele,space,iter,Lx,Ly,Lz)
open (9,file='plot',status = 'replace')
call plot_anim(n_ele,space,iter,Lx,Ly,Lz,10)
close(9)
call system('gnuplot plot')
call system('rm -rf plot')
call system('rm -rf data_frame.dat')
Expand All @@ -285,4 +329,8 @@ program Thomson

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

close(10)

call MPI_Finalize(ierr)

end program
Loading