Skip to content

Commit 00378f7

Browse files
authored
Merge pull request #307 from deepmodeling/GPU
Gpu
2 parents 5131144 + 420afda commit 00378f7

109 files changed

Lines changed: 3060339 additions & 105 deletions

Some content is hidden

Large Commits have some content hidden by default. Use the searchbox below for content that may be hidden.

.gitignore

Lines changed: 3 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -70,3 +70,6 @@ __pycache__/
7070
lib/
7171
bin/
7272
.vscode/
73+
result/
74+
*result*
75+
*profile*

Allwclean

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -25,3 +25,4 @@ wclean ./applications/solvers/dfHighSpeedFoam
2525
rm -rf src_orig/
2626
rm -rf bin/
2727
rm -rf lib/
28+
rm -rf src_gpu/build

applications/solvers/dfLowMachFoam/CMakeLists.txt

Lines changed: 10 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1,6 +1,8 @@
11
cmake_minimum_required(VERSION 3.5)
22
project(dfLowMachFoam LANGUAGES CXX)
33
FIND_PACKAGE(MPI REQUIRED)
4+
FIND_PACKAGE(OpenMP REQUIRED)
5+
FIND_PACKAGE(CUDA REQUIRED)
46

57
# Check valid thirdParty
68
if(DEFINED ENV{WM_PROJECT_DIR})
@@ -26,6 +28,8 @@ SET(SRC_ORIG $ENV{SRC_ORIG})
2628

2729
# set compilation options
2830
SET(CMAKE_EXE_LINKER_FLAGS "-fuse-ld=bfd -Xlinker --add-needed -Xlinker --no-as-needed")
31+
SET (CMAKE_C_FLAGS ${CMAKE_C_FLAGS} ${OpenMP_C_FLAGS})
32+
SET (CMAKE_CXX_FLAGS ${CMAKE_CXX_FLAGS} ${OpenMP_CXX_FLAGS})
2933

3034
SET(CMAKE_C_COMPILER g++)
3135
SET(PATH_LIB_OPENMPI "openmpi-system") # Foundation version
@@ -83,6 +87,9 @@ include_directories(
8387
${CANTERA_ROOT}/include
8488
${MPI_INCLUDE_PATH}
8589
${PROJECT_SOURCE_DIR}
90+
${CUDA_INCLUDE_DIRS}
91+
/home/runze/AmgX/AMGX/include
92+
/home/runze/deepflame-dev/src_gpu
8693
)
8794

8895
# add execution
@@ -98,6 +105,9 @@ target_link_libraries(${PROJECT_NAME}
98105
${DF_ROOT}/lib/libdfCombustionModels.so
99106
$ENV{FOAM_LIBBIN}/openmpi-system/libPstream.so
100107
${MPI_LIBRARIES}
108+
${CUDA_LIBRARIES}
109+
/home/runze/AmgX/AMGX/build/libamgxsh.so
110+
/home/runze/deepflame-dev/src_gpu/build/libdfMatrix.so
101111
)
102112

103113
if(DEFINED ENV{PYTHON_INC_DIR})

applications/solvers/dfLowMachFoam/EEqn.H

Lines changed: 116 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -1,8 +1,113 @@
11
{
22
volScalarField& he = thermo.he();
3+
#ifdef GPUSolver_
4+
start1 = std::clock();
5+
UEqn_GPU.updatePsi(&U[0][0]);
6+
UEqn_GPU.correctBoundaryConditions();
7+
U.correctBoundaryConditions();
8+
K = 0.5*magSqr(U);
9+
end1 = std::clock();
10+
time_monitor_UEqn += double(end1 - start1) / double(CLOCKS_PER_SEC);
11+
time_monitor_UEqn_correctBC += double(end1 - start1) / double(CLOCKS_PER_SEC);
312

13+
// prepare data on CPU
14+
start1 = std::clock();
15+
start2 = std::clock();
16+
// const tmp<volScalarField> alphaEff_tmp(thermo.alpha());
17+
// const volScalarField& alphaEff = alphaEff_tmp();
18+
double *alphaEff = nullptr; // tmp
19+
end2 = std::clock();
20+
int eeqn_offset = 0;
21+
int patchNum = 0;
22+
23+
forAll(he.boundaryField(), patchi)
24+
{
25+
patchNum++;
26+
const fvsPatchScalarField& pw = mesh.surfaceInterpolation::weights().boundaryField()[patchi];
27+
int patchSize = pw.size();
28+
29+
// construct gradient manually
30+
const fvPatchScalarField& hew = he.boundaryField()[patchi];
31+
const basicThermo& bThermo = basicThermo::lookupThermo(hew);
32+
const scalarField& ppw = bThermo.p().boundaryField()[patchi];
33+
fvPatchScalarField& Tw =
34+
const_cast<fvPatchScalarField&>(bThermo.T().boundaryField()[patchi]);
35+
scalarField& Tw_v = Tw;
36+
37+
Tw.evaluate();
38+
const scalarField& patchDeltaCoeff = mesh.boundary()[patchi].deltaCoeffs();
39+
const scalarField heInternal = bThermo.he(ppw, Tw, patchi)();
40+
const scalarField heBoundary = bThermo.he(ppw, Tw, mesh.boundary()[patchi].faceCells())();
41+
const scalarField patchGradMau = patchDeltaCoeff * (heInternal - heBoundary);
42+
43+
const scalarField& patchK = K.boundaryField()[patchi];
44+
// const scalarField& patchAlphaEff = alphaEff.boundaryField()[patchi]; // not H2Dcopy when use UnityLewis
45+
// const scalarField& patchGrad = he.boundaryField()[patchi].gradientBoundaryCoeffs(); // gradient_
46+
47+
// const DimensionedField<scalar, volMesh>& patchHa_ = he.boundaryField()[patchi];
48+
// const gradientEnergyFvPatchScalarField patchHa(mesh.boundary()[patchi], patchHa_);
49+
// const scalarField& patchGrad = patchHa.gradient(); // gradient_
50+
memcpy(boundary_K + eeqn_offset, &patchK[0], patchSize*sizeof(double));
51+
// memcpy(boundary_alphaEff + eeqn_offset, &patchAlphaEff[0], patchSize*sizeof(double));
52+
memcpy(boundary_gradient + eeqn_offset, &patchGradMau[0], patchSize*sizeof(double));
53+
54+
eeqn_offset += patchSize;
55+
}
56+
end1 = std::clock();
57+
time_monitor_EEqn += double(end1 - start1) / double(CLOCKS_PER_SEC);
58+
time_monitor_EEqn_mtxAssembly += double(end1 - start1) / double(CLOCKS_PER_SEC);
59+
time_monitor_EEqn_mtxAssembly_CPU_prepare += double(end1 - start1) / double(CLOCKS_PER_SEC);
60+
fprintf(stderr, "time_monitor_EEqn_mtxAssembly_CPU_prepare: %lf, build alphaEff time: %lf, patchNum: %d\n",
61+
time_monitor_EEqn_mtxAssembly_CPU_prepare,
62+
double(end2 - start2) / double(CLOCKS_PER_SEC), patchNum);
63+
64+
// prepare data on GPU
65+
start1 = std::clock();
66+
he.oldTime();
67+
K.oldTime();
68+
EEqn_GPU.prepare_data(&he.oldTime()[0], &K[0], &K.oldTime()[0], alphaEff,
69+
&dpdt[0], boundary_K, boundary_alphaEff, boundary_gradient);
70+
EEqn_GPU.sync();
71+
end1 = std::clock();
72+
time_monitor_EEqn += double(end1 - start1) / double(CLOCKS_PER_SEC);
73+
time_monitor_EEqn_mtxAssembly += double(end1 - start1) / double(CLOCKS_PER_SEC);
74+
time_monitor_EEqn_mtxAssembly_GPU_prepare += double(end1 - start1) / double(CLOCKS_PER_SEC);
75+
76+
start1 = std::clock();
77+
EEqn_GPU.initializeTimeStep();
78+
EEqn_GPU.fvm_ddt();
79+
EEqn_GPU.fvm_div();
80+
EEqn_GPU.fvm_laplacian();
81+
EEqn_GPU.fvc_ddt();
82+
EEqn_GPU.fvc_div_phi_scalar();
83+
EEqn_GPU.fvc_div_vector();
84+
EEqn_GPU.add_to_source();
85+
EEqn_GPU.sync();
86+
end1 = std::clock();
87+
time_monitor_EEqn += double(end1 - start1) / double(CLOCKS_PER_SEC);
88+
time_monitor_EEqn_mtxAssembly += double(end1 - start1) / double(CLOCKS_PER_SEC);
89+
time_monitor_EEqn_mtxAssembly_GPU_run += double(end1 - start1) / double(CLOCKS_PER_SEC);
90+
91+
// check value of mtxAssembly, no time monitor
92+
// EEqn_GPU.checkValue(true);
93+
94+
start1 = std::clock();
95+
EEqn_GPU.solve();
96+
end1 = std::clock();
97+
time_monitor_EEqn += double(end1 - start1) / double(CLOCKS_PER_SEC);
98+
time_monitor_EEqn_solve += double(end1 - start1) / double(CLOCKS_PER_SEC);
99+
100+
start1 = std::clock();
101+
EEqn_GPU.updatePsi(&he[0]);
102+
he.correctBoundaryConditions();
103+
he.write();
104+
end1 = std::clock();
105+
time_monitor_EEqn += double(end1 - start1) / double(CLOCKS_PER_SEC);
106+
time_monitor_EEqn_correctBC += double(end1 - start1) / double(CLOCKS_PER_SEC);
107+
#else
108+
start1 = std::clock();
4109
fvScalarMatrix EEqn
5-
(
110+
(
6111

7112
fvm::ddt(rho, he) + mvConvection->fvmDiv(phi, he)
8113
+ fvc::ddt(rho, K) + fvc::div(phi, K)
@@ -22,8 +127,15 @@
22127
)
23128
)
24129
);
130+
end1 = std::clock();
131+
time_monitor_EEqn += double(end1 - start1) / double(CLOCKS_PER_SEC);
132+
time_monitor_EEqn_mtxAssembly += double(end1 - start1) / double(CLOCKS_PER_SEC);
25133

26-
EEqn.relax();
27-
28-
EEqn.solve("ha");
134+
EEqn.relax();
135+
start1 = std::clock();
136+
EEqn.solve("ha");
137+
end1 = std::clock();
138+
time_monitor_EEqn += double(end1 - start1) / double(CLOCKS_PER_SEC);
139+
time_monitor_EEqn_solve += double(end1 - start1) / double(CLOCKS_PER_SEC);
140+
#endif
29141
}

applications/solvers/dfLowMachFoam/Make/options

Lines changed: 12 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -1,12 +1,15 @@
11
-include $(GENERAL_RULES)/mplibType
22

33
EXE_INC = -std=c++14 \
4+
-g \
5+
-fopenmp \
46
-Wno-unused-variable \
57
-Wno-unused-but-set-variable \
68
-Wno-old-style-cast \
79
$(PFLAGS) $(PINC) \
810
$(if $(LIBTORCH_ROOT),-DUSE_LIBTORCH,) \
911
$(if $(PYTHON_INC_DIR),-DUSE_PYTORCH,) \
12+
$(if $(AMGX_DIR),-DGPUSolver_,) \
1013
-I$(LIB_SRC)/transportModels/compressible/lnInclude \
1114
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
1215
-I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
@@ -23,7 +26,10 @@ EXE_INC = -std=c++14 \
2326
-I$(CANTERA_ROOT)/include \
2427
$(if $(LIBTORCH_ROOT),-I$(LIBTORCH_ROOT)/include,) \
2528
$(if $(LIBTORCH_ROOT),-I$(LIBTORCH_ROOT)/include/torch/csrc/api/include,) \
26-
$(PYTHON_INC_DIR)
29+
$(PYTHON_INC_DIR) \
30+
$(if $(AMGX_DIR), -I$(DF_ROOT)/src_gpu,) \
31+
$(if $(AMGX_DIR), -I/usr/local/cuda-11.6/include,) \
32+
$(if $(AMGX_DIR), -I$(AMGX_DIR)/include,)
2733

2834
EXE_LIBS = \
2935
-lcompressibleTransportModels \
@@ -44,4 +50,8 @@ EXE_LIBS = \
4450
$(if $(LIBTORCH_ROOT),-lpthread,) \
4551
$(if $(LIBTORCH_ROOT),$(DF_SRC)/dfChemistryModel/DNNInferencer/build/libDNNInferencer.so,) \
4652
$(if $(PYTHON_LIB_DIR),-L$(PYTHON_LIB_DIR),) \
47-
$(if $(PYTHON_LIB_DIR),-lpython3.8,)
53+
$(if $(PYTHON_LIB_DIR),-lpython3.8,) \
54+
$(if $(AMGX_DIR), /usr/local/cuda-11.6/lib64/libcudart.so,) \
55+
$(if $(AMGX_DIR), $(DF_ROOT)/src_gpu/build/libdfMatrix.so,) \
56+
$(if $(AMGX_DIR), $(AMGX_DIR)/build/libamgxsh.so,)
57+

applications/solvers/dfLowMachFoam/UEqn.H

Lines changed: 127 additions & 12 deletions
Original file line numberDiff line numberDiff line change
@@ -1,17 +1,132 @@
11
// Solve the Momentum equation
2+
#ifdef GPUSolver_
3+
start1 = std::clock();
4+
int offset = 0;
5+
const tmp<volScalarField> nuEff_tmp(turbulence->nuEff());
6+
const volScalarField& nuEff = nuEff_tmp();
7+
forAll(U.boundaryField(), patchi)
8+
{
9+
const scalarField& patchP = p.boundaryField()[patchi];
10+
const vectorField& patchU = U.boundaryField()[patchi];
11+
const scalarField& patchRho = rho.boundaryField()[patchi];
12+
const scalarField& patchNuEff = nuEff.boundaryField()[patchi];
213

3-
tmp<fvVectorMatrix> tUEqn
4-
(
5-
fvm::ddt(rho, U) + fvm::div(phi, U)
6-
+ turbulence->divDevRhoReff(U)
7-
);
8-
fvVectorMatrix& UEqn = tUEqn.ref();
14+
int patchSize = patchP.size();
915

10-
UEqn.relax();
16+
// boundary pressure
17+
memcpy(boundary_pressure_init+offset, &patchP[0], patchSize*sizeof(double));
18+
// boundary velocity
19+
memcpy(boundary_velocity_init+3*offset, &patchU[0][0], 3*patchSize*sizeof(double));
20+
// boundary nuEff
21+
memcpy(boundary_nuEff_init+offset, &patchNuEff[0], patchSize*sizeof(double));
22+
// boundary rho
23+
memcpy(boundary_rho_init+offset, &patchRho[0], patchSize*sizeof(double));
24+
offset += patchSize;
25+
}
26+
end1 = std::clock();
27+
time_monitor_UEqn += double(end1 - start1) / double(CLOCKS_PER_SEC);
28+
time_monitor_UEqn_mtxAssembly += double(end1 - start1) / double(CLOCKS_PER_SEC);
29+
time_monitor_UEqn_mtxAssembly_CPU_prepare += double(end1 - start1) / double(CLOCKS_PER_SEC);
1130

12-
if (pimple.momentumPredictor())
13-
{
14-
solve(UEqn == -fvc::grad(p));
31+
start1 = std::clock();
32+
UEqn_GPU.initializeTimeStep();
33+
U.oldTime();
34+
UEqn_GPU.fvm_ddt(&U.oldTime()[0][0]);
35+
UEqn_GPU.fvm_div(boundary_pressure_init, boundary_velocity_init, boundary_nuEff_init, boundary_rho_init);
36+
UEqn_GPU.fvc_grad(&p[0]);
37+
UEqn_GPU.fvc_grad_vector();
38+
UEqn_GPU.dev2T();
39+
UEqn_GPU.fvc_div_tensor(&nuEff[0]);
40+
UEqn_GPU.fvm_laplacian();
41+
UEqn_GPU.sync();
42+
end1 = std::clock();
43+
time_monitor_UEqn += double(end1 - start1) / double(CLOCKS_PER_SEC);
44+
time_monitor_UEqn_mtxAssembly += double(end1 - start1) / double(CLOCKS_PER_SEC);
45+
time_monitor_UEqn_mtxAssembly_GPU_run += double(end1 - start1) / double(CLOCKS_PER_SEC);
46+
47+
// start2 = std::clock();
48+
// fvVectorMatrix turb_source
49+
// (
50+
// turbulence->divDevRhoReff(U)
51+
// );
52+
// end2 = std::clock();
53+
// time_monitor_CPU += double(end2 - start2) / double(CLOCKS_PER_SEC);
54+
55+
// UEqn_GPU.add_fvMatrix(&turb_source.lower()[0], &turb_source.diag()[0], &turb_source.upper()[0], &turb_source.source()[0][0]);
56+
// end1 = std::clock();
57+
// time_monitor_UEqn += double(end1 - start1) / double(CLOCKS_PER_SEC);
58+
// time_monitor_UEqn_mtxAssembly += double(end1 - start1) / double(CLOCKS_PER_SEC);
59+
60+
// check value
61+
// U.oldTime();
62+
// tmp<fvVectorMatrix> tUEqn
63+
// (
64+
// fvm::ddt(rho, U)
65+
// +
66+
// fvm::div(phi, U)
67+
// +
68+
// turbulence->divDevRhoReff(U)
69+
// == -fvc::grad(p)
70+
// );
71+
// fvVectorMatrix& UEqn = tUEqn.ref();
72+
// printf("b_cpu = %e\n", UEqn.source()[1][1]);
73+
// forAll(U.boundaryField(), patchi){
74+
// labelUList sub_boundary = mesh.boundary()[patchi].faceCells();
75+
// forAll(sub_boundary, i){
76+
// if (sub_boundary[i] == 1){
77+
// printf("b_cpu_bou = %e\n", UEqn.boundaryCoeffs()[patchi][i][1]);
78+
// printf("patchi = %d, i = %d\n", patchi, i);
79+
// }
80+
// }
81+
// }
82+
// if (pimple.momentumPredictor())
83+
// {
84+
// solve(UEqn);
85+
// Info << "U_CPU\n" << U << endl;
86+
// K = 0.5*magSqr(U);
87+
// }
88+
// UEqn_GPU.checkValue(true);
89+
#else
90+
start1 = std::clock();
91+
tmp<fvVectorMatrix> tUEqn
92+
(
93+
fvm::ddt(rho, U) + fvm::div(phi, U)
94+
+ turbulence->divDevRhoReff(U)
95+
== -fvc::grad(p)
96+
);
97+
fvVectorMatrix& UEqn = tUEqn.ref();
98+
99+
end1 = std::clock();
100+
time_monitor_UEqn += double(end1 - start1) / double(CLOCKS_PER_SEC);
101+
time_monitor_UEqn_mtxAssembly += double(end1 - start1) / double(CLOCKS_PER_SEC);
102+
103+
UEqn.relax();
104+
start1 = std::clock();
105+
if (pimple.momentumPredictor())
106+
{
107+
solve(UEqn);
108+
109+
K = 0.5*magSqr(U);
110+
}
111+
end1 = std::clock();
112+
time_monitor_UEqn += double(end1 - start1) / double(CLOCKS_PER_SEC);
113+
time_monitor_UEqn_solve += double(end1 - start1) / double(CLOCKS_PER_SEC);
114+
#endif
115+
116+
// start1 = std::clock();
117+
// // // std::thread t(&dfMatrix::solve, &UEqn_GPU);
118+
// UEqn_GPU.solve();
119+
// end1 = std::clock();
120+
// time_monitor_UEqn += double(end1 - start1) / double(CLOCKS_PER_SEC);
121+
// time_monitor_UEqn_solve += double(end1 - start1) / double(CLOCKS_PER_SEC);
122+
123+
// start1 = std::clock();
124+
// // // t.join();
125+
// // UEqn_GPU.updatePsi(&U[0][0]);
126+
// K = 0.5*magSqr(U);
127+
// end1 = std::clock();
128+
// time_monitor_UEqn += double(end1 - start1) / double(CLOCKS_PER_SEC);
129+
// time_monitor_UEqn_mtxAssembly += double(end1 - start1) / double(CLOCKS_PER_SEC);
130+
// time_monitor_CPU += double(end1 - start1) / double(CLOCKS_PER_SEC);
131+
// // Info << "U_amgx = " << U << endl;
15132

16-
K = 0.5*magSqr(U);
17-
}

0 commit comments

Comments
 (0)