I am going to buy a new powerful desktop that I will use to run finite element simulations using AceFEM. I found out that the new AMD Ryzen 9 3950X obtained more points in the CPU Passmark benchmark compared to the Intel processors 1. I would like to know how well it works with AceFEM that uses linear solvers from the MKL library (both direct Pardiso and iterative CG) and performs well with Intel processors. On the other hand, the AMD processors are not as good as Intel processors in working with MKL library. For instance, it has been observed that in Matlab, which uses MKL library, the new AMD processors were not as good as Intel processors. However, after some tweaking, the performance of AMD processors increased significantly 2.
I would like to know if anyone has experience in working with AceFEM with AMD processors? In particular, the new AMD processors, e.g. Ryzen 9 3900, 3900X or 3950X. Can AMD processors compete with Intel processors?
Just in case someone has AMD processors, I would like to do a simple test for which I will provide the related code.
1 https://www.cpubenchmark.net/cpu.php?cpu=AMD+Ryzen+9+3950X&id=3598
2 https://www.pugetsystems.com/labs/hpc/How-To-Use-MKL-with-AMD-Ryzen-and-Threadripper-CPU-s-Effectively-for-Python-Numpy-And-Other-Applications-1637
Edit
I have prepared a sample simulation to test the performance of the CPUs. The compression of a hyperelastic material (with Neo-Hookean elastic strain energy) is considered, where the Neumann boundary conditions are prescribed to induce the compression within 10 time-steps. It is important to test the performance with both direct and iterative solvers. The material element (AceGen code) and the simulations (for both direct and iterative solvers) are provided below.
AceGen code
<< AceGen`;
nNodes = 8;
nhdata = 9;
SMSInitialize["HEISONEO", "Environment" -> "AceFEM"];
SMSTemplate["SMSTopology" -> "H1", "SMSNoNodes" -> nNodes,
"SMSDOFGlobal" -> Table[3, nNodes],
"SMSNodeID" -> Table["D", nNodes], "SMSSymmetricTangent" -> True,
"SMSDefaultIntegrationCode" -> 11,
"SMSNoElementData" -> nhdata es$$["id", "NoIntPoints"],
"SMSDomainDataNames" -> {"Ee -elastic modulus",
"ν -Poisson ratio"}];
SMSStandardModule["Tangent and residual"];
initialization1[] := (
Xi ⊢ SMSReal[Table[nd$$[i, "X", j], {i, nNodes}, {j, 3}]];
ui ⊢
SMSReal[Table[nd$$[i, "at", j], {i, nNodes}, {j, 3}]];
{Ee, ν} ⊢
SMSReal[Table[es$$["Data", i], {i, Length[SMSDomainDataNames]}]];
);
initialization1[];
SMSDo[Ig, 1, SMSInteger[es$$["id", "NoIntPoints"]]];
initialization2[] := (
Ξ = {ξ, η, ζ} ⊢
Table[SMSReal[es$$["IntPoints", i, Ig]], {i, 3}];
Nodeξηζ = {{-1 , -1 , -1} , {1 , -1 , -1} , {1 ,
1 , -1} , {-1 , 1 , -1} , {-1 , -1 , 1} , {1 , -1 , 1} , {1 ,
1 , 1} , {-1 , 1 , 1}};
Ni ⊨
Table[1/8 (1 + ξ Nodeξηζ[[i ,
1]]) (1 + η Nodeξηζ[[i,
2]]) (1 + ζ Nodeξηζ[[i, 3]]) , {i ,
1 , 8}];
X ⊨ SMSFreeze[Ni.Xi];
u ⊨ Ni.ui;
Jg ⊨ SMSD[X, Ξ];
Jgd ⊨ Det[Jg];
\[DoubleStruckCapitalH] ⊨
SMSD[u, X, "Dependency" -> {Ξ, X, SMSInverse[Jg]}];
Ii ⊨ IdentityMatrix[3];
F ⊨ SMSFreeze[Ii + \[DoubleStruckCapitalH]];
SMSFreeze[Fe, F, "Ignore" -> PossibleZeroQ];
SMSExport[Flatten[Fe], ed$$["Data", (Ig - 1) nhdata + #] &];
Ce ⊨ Transpose[Fe].Fe;
be ⊨ Fe.Transpose[Fe];
Je ⊨ SMSSqrt[Det[be]];
{μ, κ} ⊨ SMSHookeToBulk[Ee, ν];
W ⊨ κ/2 (1/2 (Je^2 - 1) - Log[Je]) + μ/
2 (Tr[Je^(-2/3) be] - 3);
wgp ⊨ SMSReal[es$$["IntPoints", 4, Ig]];
);
initialization2[];
pe = Flatten[ui];
SMSDo[
Rg ⊨ Jgd wgp (SMSD[W, pe, i]);
SMSExport[SMSResidualSign Rg, p$$[i], "AddIn" -> True];
SMSDo[
Kg ⊨ SMSD[Rg, pe, j];
SMSExport[Kg, s$$[i, j], "AddIn" -> True];
, {j, SMSNoDOFGlobal}];
, {i, SMSNoDOFGlobal}];
SMSEndDo[];
SMSStandardModule["Postprocessing"];
initialization1[];
SMSNPostNames = {"DeformedMeshX", "DeformedMeshY", "DeformedMeshZ"};
SMSExport[Table[ui[[i]], {i, 8}], npost$$];
SMSWrite[];
SMTMakeDll[];
Test 1: Direct solver
<< AceFEM`;
width = 200;
depth = 200;
height = 200;
nx = 50;
ny = 50;
nz = 50;
points = {{0, 0, 0}, {width, 0, 0}, {width, depth, 0}, {0, depth,
0}, {0, 0, height}, {width, 0, height}, {width, depth,
height}, {0, depth, height}};
SMTInputData[];
SMTAddDomain[{"A", "HEISONEO", {"Ee *" -> 107.2, "ν *" -> 0.36}}];
SMTAddMesh[Hexahedron[points], "A", "H1", {nx, ny, nz}];
SMTAddEssentialBoundary[
Polygon[{{0, 0, 0}, {width, 0, 0}, {width, depth, 0}, {0, depth,
0}}, "D"], 3 -> 0];
SMTAddEssentialBoundary[
Polygon[{{0, 0, 0}, {0, depth, 0}, {0, depth, height}, {0, 0,
height}}, "D"], 1 -> 0];
SMTAddEssentialBoundary[
Polygon[{{0, 0, 0}, {width, 0, 0}, {width, 0, height}, {0, 0,
height}}, "D"], 2 -> 0];
SMTAddNaturalBoundary[
Polygon[{{0, 0, height}, {width/2, 0, height}, {width/2, depth/2,
height}, {0, depth/2, height}}, "D"], 3 -> -18];
SMTAnalysis["Output" -> "testSolver5Mat2.out"];
SMTSetSolver[5, "MatrixType" -> 2];
velocity = 10;
fd = {{0, 0}};
λf[t_] := velocity t;
SMTNextStep["Δt" -> 4/velocity,
"λ[t]" -> λf];
While[
While[
step =
SMTConvergence[1*10^-12,
16, {"Adaptive Time", 8, 4/velocity, 4/velocity, 40/velocity}],
SMTNewtonIteration[];];
If[step[[4]] === "MinBound",
Print["Error: Δt < Δtmin"]];
If[step[[4]] === "MinBound",
SMTStatusReport[
"ΔT<\!\(\*SubscriptBox[\(ΔT\), \(min\
\)]\)"];];
step[[3]], If[step[[1]], SMTStepBack[];];
SMTNextStep["Δt" -> step[[2]],
"λ[t]" -> λf];
];
SMTSimulationReport[]
Test2: Iterative solver
<< AceFEM`;
width = 200;
depth = 200;
height = 200;
nx = 50;
ny = 50;
nz = 50;
points = {{0, 0, 0}, {width, 0, 0}, {width, depth, 0}, {0, depth,
0}, {0, 0, height}, {width, 0, height}, {width, depth,
height}, {0, depth, height}};
SMTInputData[];
SMTAddDomain[{"A", "HEISONEO", {"Ee *" -> 107.2, "ν *" -> 0.36}}];
SMTAddMesh[Hexahedron[points], "A", "H1", {nx, ny, nz}];
SMTAddEssentialBoundary[
Polygon[{{0, 0, 0}, {width, 0, 0}, {width, depth, 0}, {0, depth,
0}}, "D"], 3 -> 0];
SMTAddEssentialBoundary[
Polygon[{{0, 0, 0}, {0, depth, 0}, {0, depth, height}, {0, 0,
height}}, "D"], 1 -> 0];
SMTAddEssentialBoundary[
Polygon[{{0, 0, 0}, {width, 0, 0}, {width, 0, height}, {0, 0,
height}}, "D"], 2 -> 0];
SMTAddNaturalBoundary[
Polygon[{{0, 0, height}, {width/2, 0, height}, {width/2, depth/2,
height}, {0, depth/2, height}}, "D"], 3 -> -18];
SMTAnalysis["Output" -> "testSolver6Mat2.out"];
SMTSetSolver[6, "MatrixType" -> 2, "IterativeSolverType" -> 2,
"Preconditioner" -> 3];
velocity = 10;
fd = {{0, 0}};
λf[t_] := velocity t;
SMTNextStep["Δt" -> 4/velocity,
"λ[t]" -> λf];
While[
While[
step =
SMTConvergence[1*10^-12,
16, {"Adaptive Time", 8, 4/velocity, 4/velocity, 40/velocity}],
SMTNewtonIteration[];];
If[step[[4]] === "MinBound",
Print["Error: Δt < Δtmin"]];
If[step[[4]] === "MinBound",
SMTStatusReport[
"ΔT<\!\(\*SubscriptBox[\(ΔT\), \(min\
\)]\)"];];
step[[3]], If[step[[1]], SMTStepBack[];];
SMTNextStep["Δt" -> step[[2]],
"λ[t]" -> λf];
];
SMTSimulationReport[]
I have performed the simulations with an Intel(R) Core(TM) i7-6950X CPU, Windows 64 bit, AceFEM version 6.823, and the simulation reports for direct and iterative solver are as follows.
Update (2/1/2023)
We recently bought a new desktop equipped with Intel Core i7-13700k CPU. It has 8 performance cores (16 threads) and 8 efficient cores (8 threads). The results of the test are as follows:
Direct solver:
Direct solver (affinity tweaked): Efficient cores not involved
Iterative solver:
Iterative solver (affinity tweaked): Efficient cores not involved
It shows that when only performance cores are at play, the computational efficiency (mostly of linear solver) exhibits 10%-20% improvement.
