From 4160c4fdb45e0be4c131a090b0ca8003cc9d981b Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sun, 15 Aug 2021 08:04:27 +0200 Subject: [PATCH] fix for parallel HDF5 if filters are applied, writing from one process does not work if the file is opened for parallel write --- PRIVATE | 2 +- cmake/Compiler-GNU.cmake | 4 ++-- src/DAMASK_interface.f90 | 4 ++-- src/grid/grid_mech_FEM.f90 | 22 ++++++++++-------- src/grid/grid_mech_spectral_basic.f90 | 24 ++++++++++++-------- src/grid/grid_mech_spectral_polarisation.f90 | 22 ++++++++++-------- 6 files changed, 45 insertions(+), 33 deletions(-) diff --git a/PRIVATE b/PRIVATE index ad4a685d4..7d783328f 160000 --- a/PRIVATE +++ b/PRIVATE @@ -1 +1 @@ -Subproject commit ad4a685d454271a400914c0334c017e3ac0dfc99 +Subproject commit 7d783328ff5deac313eb7de951d51d19b5883b84 diff --git a/cmake/Compiler-GNU.cmake b/cmake/Compiler-GNU.cmake index b9a7406dc..d0e7d81a2 100644 --- a/cmake/Compiler-GNU.cmake +++ b/cmake/Compiler-GNU.cmake @@ -12,9 +12,9 @@ endif () if (OPTIMIZATION STREQUAL "OFF") set (OPTIMIZATION_FLAGS "-O0") elseif (OPTIMIZATION STREQUAL "DEFENSIVE") - set (OPTIMIZATION_FLAGS "-O2") + set (OPTIMIZATION_FLAGS "-O2 -mtune=generic") elseif (OPTIMIZATION STREQUAL "AGGRESSIVE") - set (OPTIMIZATION_FLAGS "-O3 -ffast-math -funroll-loops -ftree-vectorize") + set (OPTIMIZATION_FLAGS "-O3 -march=native -ffast-math -funroll-loops -ftree-vectorize") endif () set (STANDARD_CHECK "-std=f2018 -pedantic-errors" ) diff --git a/src/DAMASK_interface.f90 b/src/DAMASK_interface.f90 index ca3179afc..7bfe93f9a 100644 --- a/src/DAMASK_interface.f90 +++ b/src/DAMASK_interface.f90 @@ -189,10 +189,10 @@ subroutine DAMASK_interface_init if (len_trim(workingDirArg) > 0) & print'(a)', ' Working dir argument: '//trim(workingDirArg) print'(a)', ' Geometry argument: '//trim(geometryArg) - print'(a)', ' Loadcase argument: '//trim(loadcaseArg) + print'(a)', ' Load case argument: '//trim(loadcaseArg) print'(/,a)', ' Working directory: '//getCWD() print'(a)', ' Geometry file: '//interface_geomFile - print'(a)', ' Loadcase file: '//interface_loadFile + print'(a)', ' Load case file: '//interface_loadFile print'(a)', ' Solver job name: '//getSolverJobName() if (interface_restartInc > 0) & print'(a,i6.6)', ' Restart from increment: ', interface_restartInc diff --git a/src/grid/grid_mech_FEM.f90 b/src/grid/grid_mech_FEM.f90 index 9cccd2fc0..ca616dffb 100644 --- a/src/grid/grid_mech_FEM.f90 +++ b/src/grid/grid_mech_FEM.f90 @@ -446,22 +446,26 @@ subroutine grid_mechanical_FEM_restartWrite fileHandle = HDF5_openFile(getSolverJobName()//'_restart.hdf5','w') groupHandle = HDF5_addGroup(fileHandle,'solver') - - call HDF5_write(P_aim,groupHandle,'P_aim',.false.) - call HDF5_write(F_aim,groupHandle,'F_aim',.false.) - call HDF5_write(F_aim_lastInc,groupHandle,'F_aim_lastInc',.false.) - call HDF5_write(F_aimDot,groupHandle,'F_aimDot',.false.) call HDF5_write(F,groupHandle,'F') call HDF5_write(F_lastInc,groupHandle,'F_lastInc') call HDF5_write(u_current,groupHandle,'u') call HDF5_write(u_lastInc,groupHandle,'u_lastInc') - - call HDF5_write(C_volAvg,groupHandle,'C_volAvg',.false.) - call HDF5_write(C_volAvgLastInc,groupHandle,'C_volAvgLastInc',.false.) - call HDF5_closeGroup(groupHandle) call HDF5_closeFile(fileHandle) + if (worldrank == 0) then + fileHandle = HDF5_openFile(getSolverJobName()//'_restart.hdf5','a',.false.) + groupHandle = HDF5_openGroup(fileHandle,'solver') + call HDF5_write(P_aim,groupHandle,'P_aim',.false.) + call HDF5_write(F_aim,groupHandle,'F_aim',.false.) + call HDF5_write(F_aim_lastInc,groupHandle,'F_aim_lastInc',.false.) + call HDF5_write(F_aimDot,groupHandle,'F_aimDot',.false.) + call HDF5_write(C_volAvg,groupHandle,'C_volAvg',.false.) + call HDF5_write(C_volAvgLastInc,groupHandle,'C_volAvgLastInc',.false.) + call HDF5_closeGroup(groupHandle) + call HDF5_closeFile(fileHandle) + endif + call DMDAVecRestoreArrayF90(mechanical_grid,solution_current,u_current,ierr) CHKERRQ(ierr) call DMDAVecRestoreArrayF90(mechanical_grid,solution_lastInc,u_lastInc,ierr) diff --git a/src/grid/grid_mech_spectral_basic.f90 b/src/grid/grid_mech_spectral_basic.f90 index 6f382d639..b033a7b29 100644 --- a/src/grid/grid_mech_spectral_basic.f90 +++ b/src/grid/grid_mech_spectral_basic.f90 @@ -389,21 +389,25 @@ subroutine grid_mechanical_spectral_basic_restartWrite fileHandle = HDF5_openFile(getSolverJobName()//'_restart.hdf5','w') groupHandle = HDF5_addGroup(fileHandle,'solver') - - call HDF5_write(P_aim,groupHandle,'P_aim',.false.) - call HDF5_write(F_aim,groupHandle,'F_aim',.false.) - call HDF5_write(F_aim_lastInc,groupHandle,'F_aim_lastInc',.false.) - call HDF5_write(F_aimDot,groupHandle,'F_aimDot',.false.) call HDF5_write(F,groupHandle,'F') call HDF5_write(F_lastInc,groupHandle,'F_lastInc') - - call HDF5_write(C_volAvg,groupHandle,'C_volAvg',.false.) - call HDF5_write(C_volAvgLastInc,groupHandle,'C_volAvgLastInc',.false.) - call HDF5_write(C_minMaxAvg,groupHandle,'C_minMaxAvg',.false.) - call HDF5_closeGroup(groupHandle) call HDF5_closeFile(fileHandle) + if (worldrank == 0) then + fileHandle = HDF5_openFile(getSolverJobName()//'_restart.hdf5','a',.false.) + groupHandle = HDF5_openGroup(fileHandle,'solver') + call HDF5_write(P_aim,groupHandle,'P_aim',.false.) + call HDF5_write(F_aim,groupHandle,'F_aim',.false.) + call HDF5_write(F_aim_lastInc,groupHandle,'F_aim_lastInc',.false.) + call HDF5_write(F_aimDot,groupHandle,'F_aimDot',.false.) + call HDF5_write(C_volAvg,groupHandle,'C_volAvg',.false.) + call HDF5_write(C_volAvgLastInc,groupHandle,'C_volAvgLastInc',.false.) + call HDF5_write(C_minMaxAvg,groupHandle,'C_minMaxAvg',.false.) + call HDF5_closeGroup(groupHandle) + call HDF5_closeFile(fileHandle) + endif + if (num%update_gamma) call utilities_saveReferenceStiffness call DMDAVecRestoreArrayF90(da,solution_vec,F,ierr); CHKERRQ(ierr) diff --git a/src/grid/grid_mech_spectral_polarisation.f90 b/src/grid/grid_mech_spectral_polarisation.f90 index 7a4ae7595..8258ad43d 100644 --- a/src/grid/grid_mech_spectral_polarisation.f90 +++ b/src/grid/grid_mech_spectral_polarisation.f90 @@ -445,22 +445,26 @@ subroutine grid_mechanical_spectral_polarisation_restartWrite fileHandle = HDF5_openFile(getSolverJobName()//'_restart.hdf5','w') groupHandle = HDF5_addGroup(fileHandle,'solver') - - call HDF5_write(F_aim,groupHandle,'P_aim',.false.) - call HDF5_write(F_aim,groupHandle,'F_aim',.false.) - call HDF5_write(F_aim_lastInc,groupHandle,'F_aim_lastInc',.false.) - call HDF5_write(F_aimDot,groupHandle,'F_aimDot',.false.) call HDF5_write(F,groupHandle,'F') call HDF5_write(F_lastInc,groupHandle,'F_lastInc') call HDF5_write(F_tau,groupHandle,'F_tau') call HDF5_write(F_tau_lastInc,groupHandle,'F_tau_lastInc') - - call HDF5_write(C_volAvg,groupHandle,'C_volAvg',.false.) - call HDF5_write(C_volAvgLastInc,groupHandle,'C_volAvgLastInc',.false.) - call HDF5_closeGroup(groupHandle) call HDF5_closeFile(fileHandle) + if (worldrank == 0) then + fileHandle = HDF5_openFile(getSolverJobName()//'_restart.hdf5','a',.false.) + groupHandle = HDF5_openGroup(fileHandle,'solver') + call HDF5_write(F_aim,groupHandle,'P_aim',.false.) + call HDF5_write(F_aim,groupHandle,'F_aim',.false.) + call HDF5_write(F_aim_lastInc,groupHandle,'F_aim_lastInc',.false.) + call HDF5_write(F_aimDot,groupHandle,'F_aimDot',.false.) + call HDF5_write(C_volAvg,groupHandle,'C_volAvg',.false.) + call HDF5_write(C_volAvgLastInc,groupHandle,'C_volAvgLastInc',.false.) + call HDF5_closeGroup(groupHandle) + call HDF5_closeFile(fileHandle) + endif + if(num%update_gamma) call utilities_saveReferenceStiffness call DMDAVecRestoreArrayF90(da,solution_vec,FandF_tau,ierr); CHKERRQ(ierr)