From e229fc02bf76598082ad54bab01a12421ba14f18 Mon Sep 17 00:00:00 2001 From: Young Date: Wed, 28 May 2025 21:26:13 -0600 Subject: [PATCH 01/28] rename greek variables, rename duplicate k fn, comment unused fns --- pvade/structure/ModalAnalysis.py | 50 ++++++++++++++++++-------------- 1 file changed, 28 insertions(+), 22 deletions(-) diff --git a/pvade/structure/ModalAnalysis.py b/pvade/structure/ModalAnalysis.py index 732943c..6216350 100644 --- a/pvade/structure/ModalAnalysis.py +++ b/pvade/structure/ModalAnalysis.py @@ -93,13 +93,17 @@ def __init__(self, domain, structural_analysis, params): # Define structural properties self.E = params.structure.elasticity_modulus # 1.0e9 - self.ν = params.structure.poissons_ratio # 0.3 - self.μ = self.E / (2.0 * (1.0 + self.ν)) - self.λ = self.E * self.ν / ((1.0 + self.ν) * (1.0 - 2.0 * self.ν)) + self.poissons_ratio = params.structure.poissons_ratio # 0.3 + self.lame_mu = self.E / (2.0 * (1.0 + self.poissons_ratio)) + self.lame_lambda = ( + self.E + * self.poissons_ratio + / ((1.0 + self.poissons_ratio) * (1.0 - 2.0 * self.poissons_ratio)) + ) if self.rank == 0: print( - f"mu = {self.μ} lambda = {self.λ} E = {self.E} nu = {self.ν} density = {self.rho.value}" + f"mu = {self.lame_mu} lambda = {self.lame_lambda} E = {self.E} nu = {self.poissons_ratio} density = {self.rho.value}" ) # time step @@ -334,42 +338,44 @@ def build_forms(self, domain, params): # # Stress tensor # def sigma(r): - # return dolfinx.fem.form(2.0*self.μ*ufl.sym(ufl.grad(r)) + self.λ *ufl.tr(ufl.sym(ufl.grad(r)))*ufl.Identity(len(r))) + # return dolfinx.fem.form(2.0*self.lame_mu*ufl.sym(ufl.grad(r)) + self.lame_lambda *ufl.tr(ufl.sym(ufl.grad(r)))*ufl.Identity(len(r))) # # Mass form # def m(u, u_): # return dolfinx.fem.form(self.rho*ufl.inner(u, u_)*ufl.dx) # # Elastic stiffness form - # def k(u, u_): + # def k_nominal(u, u_): # return dolfinx.fem.form(ufl.inner(sigma(u), ufl.sym(ufl.grad(u_)))*ufl.dx) # # Rayleigh damping form # def c(u, u_): - # return dolfinx.fem.form(self.eta_m*m(u, u_) + self.eta_k*k(u, u_)) + # return dolfinx.fem.form(self.eta_m*m(u, u_) + self.eta_k*k_nominal(u, u_)) # # Work of external forces # def Wext(u_): # return ufl.dot(u_, self.f)*self.ds #dss(3) def sigma(r): - return self.λ * ufl.nabla_div(r) * ufl.Identity( + return self.lame_lambda * ufl.nabla_div(r) * ufl.Identity( len(r) - ) + 2 * self.μ * ufl.sym(ufl.grad(r)) + ) + 2 * self.lame_mu * ufl.sym(ufl.grad(r)) def m(u, u_): return self.rho * ufl.inner(u, u_) * ufl.dx - def k(u, u_): + def k_cauchy(u, u_): return ufl.inner(sigma(u), ufl.grad(u_)) * ufl.dx - def k(u, u_): + def k_nominal(u, u_): # return ufl.inner(sigma(u),ufl.grad(u_))*ufl.dx # updated lagrangian form F = ufl.grad(u) + ufl.Identity(len(u)) E = 0.5 * (F.T * F - ufl.Identity(len(u))) - # S = self.λ *ufl.tr(E)*ufl.Identity(len(u)) + 2*self.μ * (E - ufl.tr(E)*ufl.Identity(len(u)) /3.0) - S = self.λ * ufl.tr(E) * ufl.Identity(len(u)) + 2 * self.μ * (E) + # S = self.lame_lambda *ufl.tr(E)*ufl.Identity(len(u)) + 2*self.lame_mu * (E - ufl.tr(E)*ufl.Identity(len(u)) /3.0) + S = self.lame_lambda * ufl.tr(E) * ufl.Identity( + len(u) + ) + 2 * self.lame_mu * (E) # return ufl.inner(F * S, ufl.grad(u_)) * ufl.dx return ufl.inner(P_(u), ufl.grad(u_)) * ufl.dx @@ -400,7 +406,7 @@ def S_(u): # return lamda * ufl.tr(E) * I + 2.0 * mu * (E - ufl.tr(E) * I / 3.0) # TODO: Why does the above form give a better result and where does it come from? - S_svk = self.λ * ufl.tr(E) * I + 2.0 * self.μ * E + S_svk = self.lame_lambda * ufl.tr(E) * I + 2.0 * self.lame_mu * E return S_svk # The first Piola–Kirchhoff stress tensor, P = F*S @@ -411,15 +417,15 @@ def P_(u): return F * S def c(u, u_): - return self.eta_m * m(u, u_) + self.eta_k * k(u, u_) + return self.eta_m * m(u, u_) + self.eta_k * k_nominal(u, u_) # self.uh_exp = dolfinx.fem.Function(self.V, name="Deformation") - def σ(v): - """Return an expression for the stress σ given a displacement field""" - return 2.0 * self.μ * ufl.sym(ufl.grad(v)) + self.λ * ufl.tr( - ufl.sym(ufl.grad(v)) - ) * ufl.Identity(len(v)) + # def σ(v): + # """Return an expression for the stress σ given a displacement field""" + # return 2.0 * self.lame_mu * ufl.sym(ufl.grad(v)) + self.lame_lambda * ufl.tr( + # ufl.sym(ufl.grad(v)) + # ) * ufl.Identity(len(v)) # source term ($f = \rho \omega^2 [x_0, \, x_1]$) # self.ω, self.ρ = 300.0, 10.0 @@ -458,7 +464,7 @@ def σ(v): self.res = ( m(self.avg(self.a_old, a_new, self.alpha_m), self.u_) + c(self.avg(self.v_old, v_new, self.alpha_f), self.u_) - + k(self.avg(self.u_old, self.u, self.alpha_f), self.u_) + + k_nominal(self.avg(self.u_old, self.u, self.alpha_f), self.u_) - self.rho * ufl.inner(self.f, self.u_) * ufl.dx - ufl.dot(ufl.dot(self.stress_predicted * J * ufl.inv(F.T), n), self.u_) * self.ds @@ -584,7 +590,7 @@ def build_nullspace(self, V): def solve(self, params, dataIO): # def σ(v): # """Return an expression for the stress σ given a displacement field""" - # return 2.0 * self.μ * ufl.sym(ufl.grad(v)) + self.λ * ufl.tr( + # return 2.0 * self.lame_mu * ufl.sym(ufl.grad(v)) + self.lame_lambda * ufl.tr( # ufl.sym(ufl.grad(v)) # ) * ufl.Identity(len(v)) From 05e288250ee070d279809dc153882c3dd61b8e9a Mon Sep 17 00:00:00 2001 From: Young Date: Wed, 28 May 2025 21:34:55 -0600 Subject: [PATCH 02/28] remove unused imports, create linear problem, solve and write tensor to structural xdmf file --- pvade/IO/DataStream.py | 55 ++++++++++++++++++++++++++++-------------- 1 file changed, 37 insertions(+), 18 deletions(-) diff --git a/pvade/IO/DataStream.py b/pvade/IO/DataStream.py index 2c158b6..d1982cd 100644 --- a/pvade/IO/DataStream.py +++ b/pvade/IO/DataStream.py @@ -1,22 +1,11 @@ -# from dolfinx import * -import numpy as np -import time -import os -import shutil -from dolfinx.io import XDMFFile, VTKFile -from mpi4py import MPI -from pathlib import Path -import pytest import dolfinx -from petsc4py import PETSc -import json +import ufl -# from dolfinx.fem import create_nonmatching_meshes_interpolation_data +import numpy as np -# hello +# from dolfinx.fem import create_nonmatching_meshes_interpolation_data -# test actions class DataStream: """Input/Output and file writing class @@ -61,7 +50,9 @@ def __init__(self, domain, flow, structure, params): f"{params.general.output_dir_sol}/solution_fluid.xdmf" ) - with XDMFFile(self.comm, self.results_filename_fluid, "w") as xdmf_file: + with dolfinx.io.XDMFFile( + self.comm, self.results_filename_fluid, "w" + ) as xdmf_file: tt = 0.0 xdmf_file.write_mesh(domain.fluid.msh) xdmf_file.write_function(flow.u_k, 0.0) @@ -78,13 +69,31 @@ def __init__(self, domain, flow, structure, params): f"{params.general.output_dir_sol}/solution_structure.xdmf" ) - with XDMFFile(self.comm, self.results_filename_structure, "w") as xdmf_file: + # Since this is the first call, build solvers that project forms on an as-saved basis + # Set up the linear problem used for the projection, cg solver and jacobi pc + petsc_options = {"ksp_type": "cg", "pc_type": "jacobi"} + self.prob_k_nominal = dolfinx.fem.petsc.LinearProblem( + ufl.lhs(structure.elasticity.k_nominal_proj), + ufl.rhs(structure.elasticity.k_nominal_proj), + bcs=[], + petsc_options=petsc_options, + ) + + # Solve for the stress tensor + sol_k_nominal = self.prob_k_nominal.solve() + structure.elasticity.internal_stress.x.array[:] = sol_k_nominal.x.array[:] + structure.elasticity.internal_stress.x.scatter_forward() + + with dolfinx.io.XDMFFile( + self.comm, self.results_filename_structure, "w" + ) as xdmf_file: tt = 0.0 xdmf_file.write_mesh(domain.structure.msh) xdmf_file.write_function(structure.elasticity.u, 0.0) xdmf_file.write_function(structure.elasticity.stress, 0.0) xdmf_file.write_function(structure.elasticity.v_old, 0.0) xdmf_file.write_function(structure.elasticity.sigma_vm_h, 0.0) + xdmf_file.write_function(structure.elasticity.internal_stress, 0.0) if self.comm.rank == 0 and self.comm.size > 1 and params.general.test: self.log_filename_structure = f"{params.general.output_dir_sol}/log_str.txt" @@ -151,7 +160,9 @@ def save_XDMF_files(self, fsi_object, domain, tt): """ if fsi_object.name == "fluid": - with XDMFFile(self.comm, self.results_filename_fluid, "a") as xdmf_file: + with dolfinx.io.XDMFFile( + self.comm, self.results_filename_fluid, "a" + ) as xdmf_file: xdmf_file.write_function(fsi_object.u_k, tt) xdmf_file.write_function(fsi_object.p_k, tt) xdmf_file.write_function(fsi_object.panel_stress, tt) @@ -160,11 +171,19 @@ def save_XDMF_files(self, fsi_object, domain, tt): xdmf_file.write_function(fsi_object.theta_k, tt) elif fsi_object.name == "structure": - with XDMFFile(self.comm, self.results_filename_structure, "a") as xdmf_file: + # Solve for the stress tensor + sol_k_nominal = self.prob_k_nominal.solve() + fsi_object.elasticity.internal_stress.x.array[:] = sol_k_nominal.x.array[:] + fsi_object.elasticity.internal_stress.x.scatter_forward() + + with dolfinx.io.XDMFFile( + self.comm, self.results_filename_structure, "a" + ) as xdmf_file: xdmf_file.write_function(fsi_object.elasticity.u, tt) xdmf_file.write_function(fsi_object.elasticity.stress, tt) xdmf_file.write_function(fsi_object.elasticity.v_old, tt) xdmf_file.write_function(fsi_object.elasticity.sigma_vm_h, tt) + xdmf_file.write_function(fsi_object.elasticity.internal_stress, tt) else: raise ValueError( From ab35b419e1bb709b684c96b38a552056669ffd21 Mon Sep 17 00:00:00 2001 From: Young Date: Wed, 28 May 2025 21:38:46 -0600 Subject: [PATCH 03/28] add form for projecting P_, remove ufl.dx from k_nominal and add directly to res, add internal stress and trial/test functions on tensor function space --- pvade/structure/ElasticityAnalysis.py | 153 +++++++------------------- 1 file changed, 38 insertions(+), 115 deletions(-) diff --git a/pvade/structure/ElasticityAnalysis.py b/pvade/structure/ElasticityAnalysis.py index 13ecba2..12c9e77 100644 --- a/pvade/structure/ElasticityAnalysis.py +++ b/pvade/structure/ElasticityAnalysis.py @@ -40,9 +40,6 @@ def __init__(self, domain, structural_analysis, params): self.rank = domain.rank self.num_procs = domain.num_procs - # # domain.structure.msh = dolfinx.mesh.refine(domain.structure.msh,None) - # # domain.structure.msh = dolfinx.mesh.refine(domain.structure.msh) - P1 = ufl.VectorElement("Lagrange", domain.structure.msh.ufl_cell(), 2) self.V = dolfinx.fem.FunctionSpace(domain.structure.msh, P1) @@ -69,85 +66,6 @@ def __init__(self, domain, structural_analysis, params): # time step self.dt_st = dolfinx.fem.Constant(domain.structure.msh, (params.structure.dt)) - # # Store the dimension of the problem for convenience - # self.ndim = domain.structure.msh.topology.dim - # self.facet_dim = self.ndim - 1 - - # # find hmin in mesh - # num_cells = domain.structure.msh.topology.index_map(self.ndim).size_local - # h = dolfinx.cpp.mesh.h(domain.structure.msh, self.ndim, range(num_cells)) - - # # This value of hmin is local to the mesh portion owned by the process - # hmin_local = np.amin(h) - - # # collect the minimum hmin from all processes - # self.hmin = np.zeros(1) - # self.comm.Allreduce(hmin_local, self.hmin, op=MPI.MIN) - # self.hmin = self.hmin[0] - - # if self.rank == 0: - # print(f"hmin on structure = {self.hmin}") - # print(f"Total num dofs on structure = {self.num_V_dofs}") - - # # Mass density - # self.rho = dolfinx.fem.Constant( - # domain.structure.msh, params.structure.rho - # ) # Constant(0.) - - # # Rayleigh damping coefficients - # self.eta_m = dolfinx.fem.Constant(domain.structure.msh, 0.0) # Constant(0.) - # self.eta_k = dolfinx.fem.Constant(domain.structure.msh, 0.0) # Constant(0.) - - # # Generalized-alpha method parameters - # self.alpha_m = dolfinx.fem.Constant(domain.structure.msh, 0.2) - # self.alpha_f = dolfinx.fem.Constant(domain.structure.msh, 0.4) - # self.gamma = 0.5 + self.alpha_f - self.alpha_m - # self.beta = (self.gamma + 0.5) ** 2 / 4.0 - - # # Define structural properties - # self.E = params.structure.elasticity_modulus # 1.0e9 - # self.ν = params.structure.poissons_ratio # 0.3 - # self.μ = self.E / (2.0 * (1.0 + self.ν)) - # self.λ = self.E * self.ν / ((1.0 + self.ν) * (1.0 - 2.0 * self.ν)) - - # if self.rank == 0: - # print( - # f"mu = {self.μ} lambda = {self.λ} E = {self.E} nu = {self.ν} density = {self.rho.value}" - # ) - - # # time step - # self.dt_st = dolfinx.fem.Constant(domain.structure.msh, (params.structure.dt)) - - # def _north_east_corner(x): - # eps = 1.0e-6 - - # # TODO: allow the probing of (x,y,z) points as something specified in the yaml file - # if params.general.geometry_module == "flag2d": - # # TEMP Hack for Turek and Hron Flag - # x1 = 0.6 - # x2 = 0.2 - # corner = [x1, x2] - # else: - # tracker_angle_rad = np.radians(params.pv_array.tracker_angle) - # x1 = 0.5 * params.pv_array.panel_chord * np.cos(tracker_angle_rad) - # x2 = 0.5 * params.pv_array.panel_thickness * np.sin(tracker_angle_rad) - # corner = [x1 - x2, 0.5 * params.pv_array.panel_span] - - # east_edge = np.logical_and(corner[0] - eps < x[0], x[0] < corner[0] + eps) - # north_edge = np.logical_and(corner[1] - eps < x[1], x[1] < corner[1] + eps) - - # north_east_corner = np.logical_and(east_edge, north_edge) - - # return north_east_corner - - # north_east_corner_facets = dolfinx.mesh.locate_entities_boundary( - # domain.structure.msh, 0, _north_east_corner - # ) - - # self.north_east_corner_dofs = dolfinx.fem.locate_dofs_topological( - # self.V, 0, north_east_corner_facets - # ) - def build_boundary_conditions(self, domain, params): """Build the boundary conditions @@ -232,6 +150,10 @@ def build_forms(self, domain, params, structure): P3 = ufl.TensorElement("Lagrange", domain.structure.msh.ufl_cell(), 2) self.T = dolfinx.fem.FunctionSpace(domain.structure.msh, P3) + self.trial_tensor = ufl.TrialFunction(self.T) + self.test_tensor = ufl.TestFunction(self.T) + self.internal_stress = dolfinx.fem.Function(self.T, name="stress_structure") + self.stress = dolfinx.fem.Function(self.T, name="stress_fluid") self.stress_old = dolfinx.fem.Function(self.T, name="stress_fluid_old") self.stress_predicted = dolfinx.fem.Function( @@ -255,47 +177,41 @@ def build_forms(self, domain, params, structure): # dss = ufl.ds(subdomain_data=boundary_subdomains) - # # Stress tensor # def sigma(r): - # return dolfinx.fem.form(2.0*self.μ*ufl.sym(ufl.grad(r)) + self.λ *ufl.tr(ufl.sym(ufl.grad(r)))*ufl.Identity(len(r))) + # return dolfinx.fem.form(2.0*self.lame_mu*ufl.sym(ufl.grad(r)) + self.lame_lambda *ufl.tr(ufl.sym(ufl.grad(r)))*ufl.Identity(len(r))) # # Mass form # def m(u, u_): # return dolfinx.fem.form(self.rho*ufl.inner(u, u_)*ufl.dx) # # Elastic stiffness form - # def k(u, u_): + # def k_nominal(u, u_): # return dolfinx.fem.form(ufl.inner(sigma(u), ufl.sym(ufl.grad(u_)))*ufl.dx) # # Rayleigh damping form # def c(u, u_): - # return dolfinx.fem.form(self.eta_m*m(u, u_) + self.eta_k*k(u, u_)) + # return dolfinx.fem.form(self.eta_m*m(u, u_) + self.eta_k*k_nominal(u, u_)) # # Work of external forces # def Wext(u_): # return ufl.dot(u_, self.f)*self.ds #dss(3) - def sigma(r): - return structure.λ * ufl.nabla_div(r) * ufl.Identity( - len(r) - ) + 2 * structure.μ * ufl.sym(ufl.grad(r)) + # def sigma(r): + # return structure.lame_lambda * ufl.nabla_div(r) * ufl.Identity( + # len(r) + # ) + 2 * structure.lame_mu * ufl.sym(ufl.grad(r)) def m(u, u_): - return structure.rho * ufl.inner(u, u_) * ufl.dx + return structure.rho * ufl.inner(u, u_) - def k(u, u_): - return ufl.inner(sigma(u), ufl.grad(u_)) * ufl.dx + def c(u, u_): + return self.eta_m * m(u, u_) + self.eta_k * k_nominal(u, u_) - def k(u, u_): - # return ufl.inner(sigma(u),ufl.grad(u_))*ufl.dx - # updated lagrangian form - F = ufl.grad(u) + ufl.Identity(len(u)) - E = 0.5 * (F.T * F - ufl.Identity(len(u))) - # S = self.λ *ufl.tr(E)*ufl.Identity(len(u)) + 2*self.μ * (E - ufl.tr(E)*ufl.Identity(len(u)) /3.0) - S = structure.λ * ufl.tr(E) * ufl.Identity(len(u)) + 2 * structure.μ * (E) + # def k_cauchy(u, u_): + # return ufl.inner(sigma(u), ufl.grad(u_)) - # return ufl.inner(F * S, ufl.grad(u_)) * ufl.dx - return ufl.inner(P_(u), ufl.grad(u_)) * ufl.dx + def k_nominal(u, u_): + return ufl.inner(P_(u), ufl.grad(u_)) # The deformation gradient, F = I + dy/dX def F_(u): @@ -323,7 +239,7 @@ def S_(u): # return lamda * ufl.tr(E) * I + 2.0 * mu * (E - ufl.tr(E) * I / 3.0) # TODO: Why does the above form give a better result and where does it come from? - S_svk = structure.λ * ufl.tr(E) * I + 2.0 * structure.μ * E + S_svk = structure.lame_lambda * ufl.tr(E) * I + 2.0 * structure.lame_mu * E return S_svk # The first Piola–Kirchhoff stress tensor, P = F*S @@ -333,16 +249,13 @@ def P_(u): # return ufl.inv(F) * S return F * S - def c(u, u_): - return self.eta_m * m(u, u_) + self.eta_k * k(u, u_) - # self.uh_exp = dolfinx.fem.Function(self.V, name="Deformation") - def σ(v): - """Return an expression for the stress σ given a displacement field""" - return 2.0 * structure.μ * ufl.sym(ufl.grad(v)) + structure.λ * ufl.tr( - ufl.sym(ufl.grad(v)) - ) * ufl.Identity(len(v)) + # def σ(v): + # """Return an expression for the stress σ given a displacement field""" + # return 2.0 * structure.lame_mu * ufl.sym(ufl.grad(v)) + structure.lame_lambda * ufl.tr( + # ufl.sym(ufl.grad(v)) + # ) * ufl.Identity(len(v)) # source term ($f = \rho \omega^2 [x_0, \, x_1]$) # self.ω, self.ρ = 300.0, 10.0 @@ -390,9 +303,9 @@ def σ(v): F = ufl.grad(self.u) + ufl.Identity(len(self.u)) J = ufl.det(F) self.res = ( - m(self.avg(self.a_old, a_new, self.alpha_m), self.u_) - + c(self.avg(self.v_old, v_new, self.alpha_f), self.u_) - + k(self.avg(self.u_old, self.u, self.alpha_f), self.u_) + m(self.avg(self.a_old, a_new, self.alpha_m), self.u_) * ufl.dx + + c(self.avg(self.v_old, v_new, self.alpha_f), self.u_) * ufl.dx + + k_nominal(self.avg(self.u_old, self.u, self.alpha_f), self.u_) * ufl.dx - structure.rho * ufl.inner(self.f, self.u_) * ufl.dx - ufl.dot(ufl.dot(self.stress_predicted * J * ufl.inv(F.T), n), self.u_) * self.ds @@ -401,6 +314,16 @@ def σ(v): # self.a = dolfinx.fem.form(ufl.lhs(res)) # self.L = dolfinx.fem.form(ufl.rhs(res)) + # Save a form to project the first Piola–Kirchhoff, P_, stress tensor in the structure + # u * v * dx = P_ * v * dx, where u and v are trial and test functions on tensor function space + F_k_nominal_proj = ufl.inner(self.trial_tensor, self.test_tensor) * ufl.dx + F_k_nominal_proj -= ( + ufl.inner(P_(self.avg(self.u_old, self.u, self.alpha_f)), self.test_tensor) + * ufl.dx + ) + + self.k_nominal_proj = F_k_nominal_proj + # self.a = dolfinx.fem.form(ufl.inner(σ(self.u), ufl.grad(self.v)) * ufl.dx) # self.L = dolfinx.fem.form( # ufl.dot(self.f, self.v) * ufl.dx @@ -518,7 +441,7 @@ def build_nullspace(self, V): def solve(self, params, dataIO, structure): # def σ(v): # """Return an expression for the stress σ given a displacement field""" - # return 2.0 * self.μ * ufl.sym(ufl.grad(v)) + self.λ * ufl.tr( + # return 2.0 * self.lame_mu * ufl.sym(ufl.grad(v)) + self.lame_lambda * ufl.tr( # ufl.sym(ufl.grad(v)) # ) * ufl.Identity(len(v)) From 03a149a2118dacdf58ed02af06e680dd20606d98 Mon Sep 17 00:00:00 2001 From: Young Date: Wed, 28 May 2025 21:39:27 -0600 Subject: [PATCH 04/28] remove unused code now moved to elasticityAnalysis, renaming greek vars --- pvade/structure/StructureMain.py | 306 +------------------------------ 1 file changed, 8 insertions(+), 298 deletions(-) diff --git a/pvade/structure/StructureMain.py b/pvade/structure/StructureMain.py index 68fd33f..5b7a1ba 100644 --- a/pvade/structure/StructureMain.py +++ b/pvade/structure/StructureMain.py @@ -42,37 +42,9 @@ def __init__(self, domain, params): self.rank = domain.rank self.num_procs = domain.num_procs - # domain.structure.msh = dolfinx.mesh.refine(domain.structure.msh,None) - # domain.structure.msh = dolfinx.mesh.refine(domain.structure.msh) - - # physics = 'ElasticityAnalysis' - - # domain_creation_module = ( - # # f"pvade.geometry.{params.general.geometry_module}.DomainCreation" - # f"pvade.structure.{physics}.Elasticity" - # ) - # try: - # # This is equivalent to "import pvade.geometry.{params.general.geometry_module}.DomainCreation as dcm" - # dcm = import_module(domain_creation_module) - # except: - # raise ValueError(f"Could not import {domain_creation_module}") - # init physics self.elasticity = Elasticity(domain, self.structural_analysis, params) - # P1 = ufl.VectorElement("Lagrange", domain.structure.msh.ufl_cell(), 2) - # self.V = dolfinx.fem.FunctionSpace(domain.structure.msh, P1) - - # self.W = dolfinx.fem.FunctionSpace( - # domain.structure.msh, ("Discontinuous Lagrange", 0) - # ) - - # self.first_call_to_solver = True - - # self.num_V_dofs = ( - # self.V.dofmap.index_map_bs * self.V.dofmap.index_map.size_global - # ) - # Store the dimension of the problem for convenience self.ndim = domain.structure.msh.topology.dim self.facet_dim = self.ndim - 1 @@ -99,25 +71,19 @@ def __init__(self, domain, params): domain.structure.msh, params.structure.rho ) # Constant(0.) - # # Rayleigh damping coefficients - # self.eta_m = dolfinx.fem.Constant(domain.structure.msh, 0.0) # Constant(0.) - # self.eta_k = dolfinx.fem.Constant(domain.structure.msh, 0.0) # Constant(0.) - - # # Generalized-alpha method parameters - # self.alpha_m = dolfinx.fem.Constant(domain.structure.msh, 0.2) - # self.alpha_f = dolfinx.fem.Constant(domain.structure.msh, 0.4) - # self.gamma = 0.5 + self.alpha_f - self.alpha_m - # self.beta = (self.gamma + 0.5) ** 2 / 4.0 - # Define structural properties self.E = params.structure.elasticity_modulus # 1.0e9 - self.ν = params.structure.poissons_ratio # 0.3 - self.μ = self.E / (2.0 * (1.0 + self.ν)) - self.λ = self.E * self.ν / ((1.0 + self.ν) * (1.0 - 2.0 * self.ν)) + self.poissons_ratio = params.structure.poissons_ratio # 0.3 + self.lame_mu = self.E / (2.0 * (1.0 + self.poissons_ratio)) + self.lame_lambda = ( + self.E + * self.poissons_ratio + / ((1.0 + self.poissons_ratio) * (1.0 - 2.0 * self.poissons_ratio)) + ) if self.rank == 0: print( - f"mu = {self.μ} lambda = {self.λ} E = {self.E} nu = {self.ν} density = {self.rho.value}" + f"mu = {self.lame_mu} lambda = {self.lame_lambda} E = {self.E} nu = {self.poissons_ratio} density = {self.rho.value}" ) def _north_east_corner(x): @@ -232,177 +198,6 @@ def build_forms(self, domain, params): """ self.elasticity.build_forms(domain, params, self) - # # Define trial and test functions for deformation - # # self.du = ufl.TrialFunction(self.V) - # self.u_ = ufl.TestFunction(self.V) - - # P3 = ufl.TensorElement("Lagrange", domain.structure.msh.ufl_cell(), 2) - # self.T = dolfinx.fem.FunctionSpace(domain.structure.msh, P3) - - # self.stress = dolfinx.fem.Function(self.T, name="stress_fluid") - # self.stress_old = dolfinx.fem.Function(self.T, name="stress_fluid_old") - # self.stress_predicted = dolfinx.fem.Function( - # self.T, name="stress_fluid_predicted" - # ) - - # self.sigma_vm_h = dolfinx.fem.Function(self.W, name="Stress") - - # # discplacement - # self.u = dolfinx.fem.Function(self.V, name="Deformation") - # self.u_old = dolfinx.fem.Function(self.V, name="Deformation_old") - # self.u_delta = dolfinx.fem.Function(self.V, name="Deformation_change") - - # # velocity - # self.v = dolfinx.fem.Function(self.V, name="Velocity") - # self.v_old = dolfinx.fem.Function(self.V, name="Velocity_old") - - # # acceleration - # self.a = dolfinx.fem.Function(self.V, name="acceleration") - # self.a_old = dolfinx.fem.Function(self.V, name="acceleration_old") - - # # dss = ufl.ds(subdomain_data=boundary_subdomains) - - # # # Stress tensor - # # def sigma(r): - # # return dolfinx.fem.form(2.0*self.μ*ufl.sym(ufl.grad(r)) + self.λ *ufl.tr(ufl.sym(ufl.grad(r)))*ufl.Identity(len(r))) - - # # # Mass form - # # def m(u, u_): - # # return dolfinx.fem.form(self.rho*ufl.inner(u, u_)*ufl.dx) - - # # # Elastic stiffness form - # # def k(u, u_): - # # return dolfinx.fem.form(ufl.inner(sigma(u), ufl.sym(ufl.grad(u_)))*ufl.dx) - - # # # Rayleigh damping form - # # def c(u, u_): - # # return dolfinx.fem.form(self.eta_m*m(u, u_) + self.eta_k*k(u, u_)) - - # # # Work of external forces - # # def Wext(u_): - # # return ufl.dot(u_, self.f)*self.ds #dss(3) - - # def sigma(r): - # return self.λ * ufl.nabla_div(r) * ufl.Identity( - # len(r) - # ) + 2 * self.μ * ufl.sym(ufl.grad(r)) - - # def m(u, u_): - # return self.rho * ufl.inner(u, u_) * ufl.dx - - # def k(u, u_): - # return ufl.inner(sigma(u), ufl.grad(u_)) * ufl.dx - - # def k(u, u_): - # # return ufl.inner(sigma(u),ufl.grad(u_))*ufl.dx - # # updated lagrangian form - # F = ufl.grad(u) + ufl.Identity(len(u)) - # E = 0.5 * (F.T * F - ufl.Identity(len(u))) - # # S = self.λ *ufl.tr(E)*ufl.Identity(len(u)) + 2*self.μ * (E - ufl.tr(E)*ufl.Identity(len(u)) /3.0) - # S = self.λ * ufl.tr(E) * ufl.Identity(len(u)) + 2 * self.μ * (E) - - # # return ufl.inner(F * S, ufl.grad(u_)) * ufl.dx - # return ufl.inner(P_(u), ufl.grad(u_)) * ufl.dx - - # # The deformation gradient, F = I + dy/dX - # def F_(u): - # I = ufl.Identity(len(u)) - # return I + ufl.grad(u) - - # # The Cauchy-Green deformation tensor, C = F.T * F - # def C_(u): - # F = F_(u) - # return F.T * F - - # # Green–Lagrange strain tensor, E = 0.5*(C - I) - # def E_(u): - # I = ufl.Identity(len(u)) - # C = C_(u) - - # return 0.5 * (C - I) - # # return 0.5 * (ufl.grad(u) + ufl.grad(u).T) - - # # The second Piola–Kirchhoff stress, S - # def S_(u): - # E = E_(u) - # I = ufl.Identity(len(u)) - - # # return lamda * ufl.tr(E) * I + 2.0 * mu * (E - ufl.tr(E) * I / 3.0) - # # TODO: Why does the above form give a better result and where does it come from? - - # S_svk = self.λ * ufl.tr(E) * I + 2.0 * self.μ * E - # return S_svk - - # # The first Piola–Kirchhoff stress tensor, P = F*S - # def P_(u): - # F = F_(u) - # S = S_(u) - # # return ufl.inv(F) * S - # return F * S - - # def c(u, u_): - # return self.eta_m * m(u, u_) + self.eta_k * k(u, u_) - - # # self.uh_exp = dolfinx.fem.Function(self.V, name="Deformation") - - # def σ(v): - # """Return an expression for the stress σ given a displacement field""" - # return 2.0 * self.μ * ufl.sym(ufl.grad(v)) + self.λ * ufl.tr( - # ufl.sym(ufl.grad(v)) - # ) * ufl.Identity(len(v)) - - # # source term ($f = \rho \omega^2 [x_0, \, x_1]$) - # # self.ω, self.ρ = 300.0, 10.0 - # # x = ufl.SpatialCoordinate(domain.structure.msh) - # # self.f = ufl.as_vector((0*self.ρ * self.ω**2 * x[0], self.ρ * self.ω**2 * x[1], 0.0)) - # # self.f_structure = dolfinx.fem.Constant( - # # domain.structure.msh, - # # (PETSc.ScalarType(0), PETSc.ScalarType(0), PETSc.ScalarType(0)), - # # ) - # # self.f = ufl.as_vector((0*self.ρ * self.ω**2 * x[0], self.ρ * self.ω**2 * x[1], 0.0)) - # # self.T = dolfinx.fem.Constant(domain.structure.msh, PETSc.ScalarType((0, 1.e-3, 0))) - # # self.f = dolfinx.fem.Constant(domain.structure.msh, PETSc.ScalarType((0,100,100))) - # self.f = dolfinx.fem.Constant( - # domain.structure.msh, - # PETSc.ScalarType( - # ( - # params.structure.body_force_x, - # params.structure.body_force_y, - # params.structure.body_force_z, - # ) - # ), - # ) - # self.ds = ufl.Measure("ds", domain=domain.structure.msh) - # n = ufl.FacetNormal(domain.structure.msh) - - # # Residual - # a_new = self.update_a( - # self.u, self.u_old, self.v_old, self.a_old, self.dt_st, self.beta, ufl=True - # ) - # v_new = self.update_v( - # a_new, self.u_old, self.v_old, self.a_old, self.dt_st, self.gamma, ufl=True - # ) - - # F = ufl.grad(self.u) + ufl.Identity(len(self.u)) - # J = ufl.det(F) - # self.res = ( - # m(self.avg(self.a_old, a_new, self.alpha_m), self.u_) - # + c(self.avg(self.v_old, v_new, self.alpha_f), self.u_) - # + k(self.avg(self.u_old, self.u, self.alpha_f), self.u_) - # - self.rho * ufl.inner(self.f, self.u_) * ufl.dx - # - ufl.dot(ufl.dot(self.stress_predicted * J * ufl.inv(F.T), n), self.u_) - # * self.ds - # ) # - Wext(self.u) - - # # self.a = dolfinx.fem.form(ufl.lhs(res)) - # # self.L = dolfinx.fem.form(ufl.rhs(res)) - - # # self.a = dolfinx.fem.form(ufl.inner(σ(self.u), ufl.grad(self.v)) * ufl.dx) - # # self.L = dolfinx.fem.form( - # # ufl.dot(self.f, self.v) * ufl.dx - # # + ufl.dot(ufl.dot(self.stress, n), self.v) * self.ds - # # ) - def _assemble_system(self, params): """Pre-assemble all LHS matrices and RHS vectors @@ -514,88 +309,3 @@ def build_nullspace(self, V): def solve(self, params, dataIO): self.elasticity.solve(params, dataIO, self) - # # def σ(v): - # # """Return an expression for the stress σ given a displacement field""" - # # return 2.0 * self.μ * ufl.sym(ufl.grad(v)) + self.λ * ufl.tr( - # # ufl.sym(ufl.grad(v)) - # # ) * ufl.Identity(len(v)) - - # if self.first_call_to_solver: - # if self.rank == 0: - # print("Starting Strutural Solution") - - # self._assemble_system(params) - - # num_its, converged = self.solver.solve(self.u) # solve the current time step - # assert converged - # self.u.x.scatter_forward() - - # # Calculate the change in the displacement (new - old) this is what moves the mesh - # self.u_delta.vector.array[:] = ( - # self.u.vector.array[:] - self.u_old.vector.array[:] - # ) - # self.u_delta.x.scatter_forward() - - # # Update old fields with new quantities - # self.update_fields( - # self.u, - # self.u_old, - # self.v_old, - # self.a_old, - # self.dt_st, - # self.beta, - # self.gamma, - # ) - - # # sigma_dev = σ(self.u) - (1 / 3) * ufl.tr(σ(self.u)) * ufl.Identity(len(self.u)) - # # sigma_vm = ufl.sqrt((3 / 2) * ufl.inner(sigma_dev, sigma_dev)) - - # # sigma_vm_expr = dolfinx.fem.Expression( - # # sigma_vm, self.W.element.interpolation_points() - # # ) - # # self.sigma_vm_h.interpolate(sigma_vm_expr) - - # # self.unorm = self.u.x.norm() - - # try: - # idx = self.north_east_corner_dofs[0] - # nw_corner_accel = self.u.vector.array[3 * idx : 3 * idx + 3].astype( - # np.float64 - # ) - # print(nw_corner_accel) - # except: - # nw_corner_accel = np.zeros(3, dtype=np.float64) - - # nw_corner_accel_global = np.zeros((self.num_procs, 3), dtype=np.float64) - - # self.comm.Gather(nw_corner_accel, nw_corner_accel_global, root=0) - - # # print(f"Acceleration at North West corner = {nw_corner_accel}") - - # if self.rank == 0: - # norm2 = np.sum(nw_corner_accel_global**2, axis=1) - # max_norm2_idx = np.argmax(norm2) - # np_accel = nw_corner_accel_global[max_norm2_idx, :] - - # accel_pos_filename = os.path.join( - # params.general.output_dir_sol, "accel_pos.csv" - # ) - - # if self.first_call_to_solver: - - # with open(accel_pos_filename, "w") as fp: - # fp.write("#x-pos,y-pos,z-pos\n") - # if self.ndim == 3: - # fp.write(f"{np_accel[0]},{np_accel[1]},{np_accel[2]}\n") - # elif self.ndim == 2: - # fp.write(f"{np_accel[0]},{np_accel[1]}\n") - - # else: - # with open(accel_pos_filename, "a") as fp: - # if self.ndim == 3: - # fp.write(f"{np_accel[0]},{np_accel[1]},{np_accel[2]}\n") - # elif self.ndim == 2: - # fp.write(f"{np_accel[0]},{np_accel[1]}\n") - - # if self.first_call_to_solver: - # self.first_call_to_solver = False From bd25235fa98659d61bd3420201317acdbdd92110 Mon Sep 17 00:00:00 2001 From: Young Date: Wed, 28 May 2025 21:44:16 -0600 Subject: [PATCH 05/28] pull updated gha file from logging branch --- .github/workflows/test_pvade.yaml | 5 ++++- 1 file changed, 4 insertions(+), 1 deletion(-) diff --git a/.github/workflows/test_pvade.yaml b/.github/workflows/test_pvade.yaml index c92cd14..3c9e05a 100644 --- a/.github/workflows/test_pvade.yaml +++ b/.github/workflows/test_pvade.yaml @@ -25,10 +25,13 @@ jobs: - name: Build Conda environment uses: conda-incubator/setup-miniconda@v3 with: + # auto-update-conda: true + python-version: "3.10" environment-file: environment.yaml + activate-environment: pvade - name: Run pytest shell: bash -l {0} - run: pytest -sv pvade/tests/ + run: PYTHONPATH=. pytest -sv pvade/tests/ # Job 2 of 2 - enforce Black formatting formatting: From 3c51786a5924ed7857f221e94c295ed778978cd1 Mon Sep 17 00:00:00 2001 From: Young Date: Thu, 29 May 2025 09:47:36 -0600 Subject: [PATCH 06/28] remove the von mises stress function --- pvade/structure/ElasticityAnalysis.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/pvade/structure/ElasticityAnalysis.py b/pvade/structure/ElasticityAnalysis.py index 12c9e77..773afee 100644 --- a/pvade/structure/ElasticityAnalysis.py +++ b/pvade/structure/ElasticityAnalysis.py @@ -160,7 +160,7 @@ def build_forms(self, domain, params, structure): self.T, name="stress_fluid_predicted" ) - self.sigma_vm_h = dolfinx.fem.Function(self.W, name="Stress") + # self.sigma_vm_h = dolfinx.fem.Function(self.W, name="Stress") # discplacement self.u = dolfinx.fem.Function(self.V, name="Deformation") From 08db451c30140c69405c0a24513ea67f7afbd44c Mon Sep 17 00:00:00 2001 From: Young Date: Thu, 29 May 2025 09:48:24 -0600 Subject: [PATCH 07/28] save the velocity and acceleration functions directly, remove saving VM stress --- pvade/IO/DataStream.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/pvade/IO/DataStream.py b/pvade/IO/DataStream.py index d1982cd..aa1c70a 100644 --- a/pvade/IO/DataStream.py +++ b/pvade/IO/DataStream.py @@ -90,9 +90,9 @@ def __init__(self, domain, flow, structure, params): tt = 0.0 xdmf_file.write_mesh(domain.structure.msh) xdmf_file.write_function(structure.elasticity.u, 0.0) + xdmf_file.write_function(structure.elasticity.v, 0.0) + xdmf_file.write_function(structure.elasticity.a, 0.0) xdmf_file.write_function(structure.elasticity.stress, 0.0) - xdmf_file.write_function(structure.elasticity.v_old, 0.0) - xdmf_file.write_function(structure.elasticity.sigma_vm_h, 0.0) xdmf_file.write_function(structure.elasticity.internal_stress, 0.0) if self.comm.rank == 0 and self.comm.size > 1 and params.general.test: @@ -180,9 +180,9 @@ def save_XDMF_files(self, fsi_object, domain, tt): self.comm, self.results_filename_structure, "a" ) as xdmf_file: xdmf_file.write_function(fsi_object.elasticity.u, tt) + xdmf_file.write_function(fsi_object.elasticity.v, tt) + xdmf_file.write_function(fsi_object.elasticity.a, tt) xdmf_file.write_function(fsi_object.elasticity.stress, tt) - xdmf_file.write_function(fsi_object.elasticity.v_old, tt) - xdmf_file.write_function(fsi_object.elasticity.sigma_vm_h, tt) xdmf_file.write_function(fsi_object.elasticity.internal_stress, tt) else: From d7134346b6c03fea722b42a9974bcb17d59c7354 Mon Sep 17 00:00:00 2001 From: Young Date: Thu, 29 May 2025 11:30:08 -0600 Subject: [PATCH 08/28] remember to save old, which after solve is equiv to new, since new is never explicitly updated, only the form is updated --- pvade/IO/DataStream.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/pvade/IO/DataStream.py b/pvade/IO/DataStream.py index aa1c70a..82134a1 100644 --- a/pvade/IO/DataStream.py +++ b/pvade/IO/DataStream.py @@ -90,8 +90,8 @@ def __init__(self, domain, flow, structure, params): tt = 0.0 xdmf_file.write_mesh(domain.structure.msh) xdmf_file.write_function(structure.elasticity.u, 0.0) - xdmf_file.write_function(structure.elasticity.v, 0.0) - xdmf_file.write_function(structure.elasticity.a, 0.0) + xdmf_file.write_function(structure.elasticity.v_old, 0.0) + xdmf_file.write_function(structure.elasticity.a_old, 0.0) xdmf_file.write_function(structure.elasticity.stress, 0.0) xdmf_file.write_function(structure.elasticity.internal_stress, 0.0) @@ -180,8 +180,8 @@ def save_XDMF_files(self, fsi_object, domain, tt): self.comm, self.results_filename_structure, "a" ) as xdmf_file: xdmf_file.write_function(fsi_object.elasticity.u, tt) - xdmf_file.write_function(fsi_object.elasticity.v, tt) - xdmf_file.write_function(fsi_object.elasticity.a, tt) + xdmf_file.write_function(fsi_object.elasticity.v_old, tt) + xdmf_file.write_function(fsi_object.elasticity.a_old, tt) xdmf_file.write_function(fsi_object.elasticity.stress, tt) xdmf_file.write_function(fsi_object.elasticity.internal_stress, tt) From 6f94f80b95f29d31fc64626379b32eb712d82c93 Mon Sep 17 00:00:00 2001 From: Young Date: Thu, 29 May 2025 11:30:51 -0600 Subject: [PATCH 09/28] name things to reduce confusion when postprocessing, remove names from unused vars --- pvade/structure/ElasticityAnalysis.py | 14 +++++++------- 1 file changed, 7 insertions(+), 7 deletions(-) diff --git a/pvade/structure/ElasticityAnalysis.py b/pvade/structure/ElasticityAnalysis.py index 773afee..4a84d98 100644 --- a/pvade/structure/ElasticityAnalysis.py +++ b/pvade/structure/ElasticityAnalysis.py @@ -163,17 +163,17 @@ def build_forms(self, domain, params, structure): # self.sigma_vm_h = dolfinx.fem.Function(self.W, name="Stress") # discplacement - self.u = dolfinx.fem.Function(self.V, name="Deformation") - self.u_old = dolfinx.fem.Function(self.V, name="Deformation_old") - self.u_delta = dolfinx.fem.Function(self.V, name="Deformation_change") + self.u = dolfinx.fem.Function(self.V, name="deformation") + self.u_old = dolfinx.fem.Function(self.V, name="deformation_old") + self.u_delta = dolfinx.fem.Function(self.V, name="deformation_change") # velocity - self.v = dolfinx.fem.Function(self.V, name="Velocity") - self.v_old = dolfinx.fem.Function(self.V, name="Velocity_old") + self.v = dolfinx.fem.Function(self.V) + self.v_old = dolfinx.fem.Function(self.V, name="velocity") # acceleration - self.a = dolfinx.fem.Function(self.V, name="acceleration") - self.a_old = dolfinx.fem.Function(self.V, name="acceleration_old") + self.a = dolfinx.fem.Function(self.V) + self.a_old = dolfinx.fem.Function(self.V, name="acceleration") # dss = ufl.ds(subdomain_data=boundary_subdomains) From 1fe88b63030295ec06d81d287074bbbb03764139 Mon Sep 17 00:00:00 2001 From: Young Date: Fri, 30 May 2025 09:43:33 -0600 Subject: [PATCH 10/28] black formatting --- pvade/IO/DataStream.py | 1 + 1 file changed, 1 insertion(+) diff --git a/pvade/IO/DataStream.py b/pvade/IO/DataStream.py index 18d3041..2b9156f 100644 --- a/pvade/IO/DataStream.py +++ b/pvade/IO/DataStream.py @@ -63,6 +63,7 @@ def flush(self): if rank == 0: print("Starting PVade Run") + class DataStream: """Input/Output and file writing class From 9c5a5deb702136f71506e1c9506cc7b51230b541 Mon Sep 17 00:00:00 2001 From: arswalid Date: Mon, 23 Jun 2025 10:42:59 -0600 Subject: [PATCH 11/28] avoid reading the mesh if analysis set to false --- pvade/geometry/MeshManager.py | 129 +++++++++++++++++----------------- 1 file changed, 65 insertions(+), 64 deletions(-) diff --git a/pvade/geometry/MeshManager.py b/pvade/geometry/MeshManager.py index be30aab..5bf6974 100644 --- a/pvade/geometry/MeshManager.py +++ b/pvade/geometry/MeshManager.py @@ -620,86 +620,87 @@ def read_mesh_files(self, read_mesh_dir, params): mesh_name = f"{sub_domain_name}_mesh.xdmf" mesh_filename = os.path.join(read_mesh_dir, mesh_name) - try: - if self.rank == 0: - print(f"Reading {sub_domain_name} mesh.") - - # Read the subdomain mesh - with dolfinx.io.XDMFFile(self.comm, mesh_filename, "r") as xdmf: - submesh = xdmf.read_mesh(name=mesh_name) - cell_tags = xdmf.read_meshtags(submesh, name="cell_tags") - ndim = submesh.topology.dim - submesh.topology.create_connectivity(ndim - 1, ndim) - facet_tags = xdmf.read_meshtags(submesh, name="facet_tags") + if params.general.fluid_analysis == True and sub_domain_name == "fluid" or params.general.structural_analysis == True and sub_domain_name == "structure": + try: + if self.rank == 0: + print(f"Reading {sub_domain_name} mesh.") + + # Read the subdomain mesh + with dolfinx.io.XDMFFile(self.comm, mesh_filename, "r") as xdmf: + submesh = xdmf.read_mesh(name=mesh_name) + cell_tags = xdmf.read_meshtags(submesh, name="cell_tags") + ndim = submesh.topology.dim + submesh.topology.create_connectivity(ndim - 1, ndim) + facet_tags = xdmf.read_meshtags(submesh, name="facet_tags") + + except: + if self.rank == 0: + print( + f"Could not find subdomain {sub_domain_name} mesh file, not reading this mesh." + ) - except: - if self.rank == 0: - print( - f"Could not find subdomain {sub_domain_name} mesh file, not reading this mesh." - ) - - else: + else: - class FSISubDomain: - pass + class FSISubDomain: + pass - if self.ndim is None: - self.ndim = submesh.topology.dim - else: - assert self.ndim == submesh.topology.dim + if self.ndim is None: + self.ndim = submesh.topology.dim + else: + assert self.ndim == submesh.topology.dim - submesh.topology.create_connectivity(self.ndim, self.ndim - 1) + submesh.topology.create_connectivity(self.ndim, self.ndim - 1) - sub_domain = FSISubDomain() + sub_domain = FSISubDomain() - sub_domain.msh = submesh - sub_domain.cell_tags = cell_tags - sub_domain.facet_tags = facet_tags + sub_domain.msh = submesh + sub_domain.cell_tags = cell_tags + sub_domain.facet_tags = facet_tags - # These elements do not need to be created when reading a mesh - # they are only used in the transfer of facet tags, and since - # those can be read directly from a file, we don't need these - sub_domain.entity_map = None - sub_domain.vertex_map = None - sub_domain.geom_map = None + # These elements do not need to be created when reading a mesh + # they are only used in the transfer of facet tags, and since + # those can be read directly from a file, we don't need these + sub_domain.entity_map = None + sub_domain.vertex_map = None + sub_domain.geom_map = None - setattr(self, sub_domain_name, sub_domain) + setattr(self, sub_domain_name, sub_domain) - if sub_domain_name == "fluid": - domain_ufl = ufl.Mesh( - self.fluid.msh.ufl_domain().ufl_coordinate_element() - ) - fluid_undeformed = dolfinx.mesh.Mesh( - self.comm, - self.fluid.msh.topology, - self.fluid.msh.geometry, - domain_ufl, - ) + if sub_domain_name == "fluid": + domain_ufl = ufl.Mesh( + self.fluid.msh.ufl_domain().ufl_coordinate_element() + ) + fluid_undeformed = dolfinx.mesh.Mesh( + self.comm, + self.fluid.msh.topology, + self.fluid.msh.geometry, + domain_ufl, + ) - fluid_undeformed.topology.create_connectivity( - self.ndim, self.ndim - 1 - ) + fluid_undeformed.topology.create_connectivity( + self.ndim, self.ndim - 1 + ) - sub_domain_undeformed = FSISubDomain() + sub_domain_undeformed = FSISubDomain() - sub_domain_undeformed.msh = fluid_undeformed - sub_domain_undeformed.cell_tags = cell_tags - sub_domain_undeformed.facet_tags = facet_tags + sub_domain_undeformed.msh = fluid_undeformed + sub_domain_undeformed.cell_tags = cell_tags + sub_domain_undeformed.facet_tags = facet_tags - # These elements do not need to be created when reading a mesh - # they are only used in the transfer of facet tags, and since - # those can be read directly from a file, we don't need these - sub_domain_undeformed.entity_map = None - sub_domain_undeformed.vertex_map = None - sub_domain_undeformed.geom_map = None + # These elements do not need to be created when reading a mesh + # they are only used in the transfer of facet tags, and since + # those can be read directly from a file, we don't need these + sub_domain_undeformed.entity_map = None + sub_domain_undeformed.vertex_map = None + sub_domain_undeformed.geom_map = None - setattr(self, "fluid_undeformed", sub_domain_undeformed) + setattr(self, "fluid_undeformed", sub_domain_undeformed) - # assert np.all(self.fluid.msh.geometry.x[:] == self.fluid_undeformed.msh.geometry.x[:]) - # assert np.shape(self.fluid.msh.geometry.x[:]) == np.shape(self.fluid_undeformed.msh.geometry.x[:]) + # assert np.all(self.fluid.msh.geometry.x[:] == self.fluid_undeformed.msh.geometry.x[:]) + # assert np.shape(self.fluid.msh.geometry.x[:]) == np.shape(self.fluid_undeformed.msh.geometry.x[:]) - if self.rank == 0: - print(f"Finished read {sub_domain_name} mesh.") + if self.rank == 0: + print(f"Finished read {sub_domain_name} mesh.") self.ndim = submesh.topology.dim From 121fc5ed7f6bb30b23e599f35563a8482d408cb9 Mon Sep 17 00:00:00 2001 From: arswalid Date: Mon, 23 Jun 2025 10:44:22 -0600 Subject: [PATCH 12/28] black formatted --- pvade/geometry/MeshManager.py | 7 ++++++- 1 file changed, 6 insertions(+), 1 deletion(-) diff --git a/pvade/geometry/MeshManager.py b/pvade/geometry/MeshManager.py index 5bf6974..27cd07e 100644 --- a/pvade/geometry/MeshManager.py +++ b/pvade/geometry/MeshManager.py @@ -620,7 +620,12 @@ def read_mesh_files(self, read_mesh_dir, params): mesh_name = f"{sub_domain_name}_mesh.xdmf" mesh_filename = os.path.join(read_mesh_dir, mesh_name) - if params.general.fluid_analysis == True and sub_domain_name == "fluid" or params.general.structural_analysis == True and sub_domain_name == "structure": + if ( + params.general.fluid_analysis == True + and sub_domain_name == "fluid" + or params.general.structural_analysis == True + and sub_domain_name == "structure" + ): try: if self.rank == 0: print(f"Reading {sub_domain_name} mesh.") From eb90544ea289c3a025347642bf2819277fe51389 Mon Sep 17 00:00:00 2001 From: arswalid Date: Tue, 24 Jun 2025 20:50:23 -0600 Subject: [PATCH 13/28] initial push to prep PR --- input/3d_cyld.yaml | 28 +++--- input/panels2d.yaml | 3 +- input/panels3d.yaml | 2 +- pvade/IO/DataStream.py | 2 +- pvade/IO/Utilities.py | 2 +- pvade/IO/input_schema.yaml | 4 +- pvade/fluid/FlowManager.py | 4 +- pvade/fluid/boundary_conditions.py | 4 +- pvade/geometry/MeshManager.py | 84 +++------------- pvade/geometry/cylinder2d/DomainCreation.py | 4 +- pvade/geometry/cylinder3d/DomainCreation.py | 4 +- pvade/structure/StructureMain.py | 19 +++- pvade/tests/conftest.py | 19 ++++ pvade/tests/test_all_inputs.py | 28 ++++++ pvade/tests/test_input_files.py | 2 +- pvade_main.py | 100 +++++++++++--------- 16 files changed, 165 insertions(+), 144 deletions(-) create mode 100644 pvade/tests/conftest.py create mode 100644 pvade/tests/test_all_inputs.py diff --git a/input/3d_cyld.yaml b/input/3d_cyld.yaml index 32ad638..8e56c2a 100644 --- a/input/3d_cyld.yaml +++ b/input/3d_cyld.yaml @@ -11,26 +11,26 @@ domain: y_max: 0.41 z_min: 0.0 z_max: 0.41 - l_char: 2.0 #.02 -pv_array: - stream_rows: 1 # not used - elevation: 1.5 # not used - stream_spacing: 7.0 # not used - panel_chord: 2.0 # not used - panel_span: 7.0 # not used - panel_thickness: 0.1 # not used - tracker_angle: 30.0 # not used + l_char: 5.0 #.02 +# pv_array: +# stream_rows: 1 # not used +# elevation: 1.5 # not used +# stream_spacing: 7.0 # not used +# panel_chord: 2.0 # not used +# panel_span: 7.0 # not used +# panel_thickness: 0.1 # not used +# tracker_angle: 30.0 # not used solver: - dt: 0.025 - t_final: 0.1 - solver1_ksp: cg + dt: 0.001 + t_final: 0.01 + solver1_ksp: gmres solver2_ksp: cg solver3_ksp: cg solver1_pc: jacobi solver2_pc: jacobi solver3_pc: jacobi - save_text_interval: 0.1 - save_xdmf_interval: 0.1 + save_text_interval: 0.001 + save_xdmf_interval: 0.001 fluid: velocity_profile_type: parabolic inflow_coeff: 16 # only necessary/used when applying a parabolic vel profile diff --git a/input/panels2d.yaml b/input/panels2d.yaml index cee63d4..7482b21 100644 --- a/input/panels2d.yaml +++ b/input/panels2d.yaml @@ -9,9 +9,10 @@ domain: x_max: 50 y_min: 0 y_max: 20 - l_char: 0.9 + l_char: 0.5 pv_array: stream_rows: 1 + span_rows: 1 elevation: 1.5 stream_spacing: 7.0 panel_chord: 2.0 diff --git a/input/panels3d.yaml b/input/panels3d.yaml index 93e4004..a4b5595 100644 --- a/input/panels3d.yaml +++ b/input/panels3d.yaml @@ -4,7 +4,7 @@ general: output_dir: output/panels3d mesh_only: False # input_mesh_dir: output/panels3d/mesh - structural_analysis: True + structural_analysis: False fluid_analysis: True domain: x_min: -20 diff --git a/pvade/IO/DataStream.py b/pvade/IO/DataStream.py index 8ebacc5..0105a00 100644 --- a/pvade/IO/DataStream.py +++ b/pvade/IO/DataStream.py @@ -74,7 +74,7 @@ class DataStream: """ - def __init__(self, domain, flow, structure, params): + def __init__(self, domain, flow, structure=None , params=None): """Initialize the DataStream object This initializes an object that manages the I/O for all of PVade. diff --git a/pvade/IO/Utilities.py b/pvade/IO/Utilities.py index 0fbb35d..250dbcf 100644 --- a/pvade/IO/Utilities.py +++ b/pvade/IO/Utilities.py @@ -23,7 +23,7 @@ def get_input_file(): return input_file_path -def write_metrics(flow, elasticity, prof_filename="profiling.txt"): +def write_metrics(flow, prof_filename="profiling.txt"): with open(prof_filename, "r") as output_file: if flow.fluid_analysis == True: # solver_line = [line for line in output_file if "(solve)" in line] diff --git a/pvade/IO/input_schema.yaml b/pvade/IO/input_schema.yaml index 0182a44..6f96806 100644 --- a/pvade/IO/input_schema.yaml +++ b/pvade/IO/input_schema.yaml @@ -118,14 +118,14 @@ properties: type: "object" properties: stream_rows: - default: 3 + default: 0 minimum: 1 maximum: 10 type: "integer" description: "The number of panel rows in the streamwise direction." units: "unitless" span_rows: - default: 1 + default: 0 minimum: 1 maximum: 10 type: "integer" diff --git a/pvade/fluid/FlowManager.py b/pvade/fluid/FlowManager.py index 5cd3f36..8b84338 100644 --- a/pvade/fluid/FlowManager.py +++ b/pvade/fluid/FlowManager.py @@ -734,7 +734,9 @@ def solve(self, domain, params, current_time): self.compute_lift_and_drag(params, current_time) - self.compute_pressure_drop_between_points(domain, params) + # Compute the pressure drop between the inlet and outlet + if params.pv_array.stream_rows > 0: + self.compute_pressure_drop_between_points(domain, params) # Update new -> old variables self.u_k2.x.array[:] = self.u_k1.x.array diff --git a/pvade/fluid/boundary_conditions.py b/pvade/fluid/boundary_conditions.py index 1d82be9..f4e6819 100644 --- a/pvade/fluid/boundary_conditions.py +++ b/pvade/fluid/boundary_conditions.py @@ -309,8 +309,8 @@ def build_velocity_boundary_conditions(domain, params, functionspace, current_ti domain, params, functionspace, current_time ) - if domain.rank == 0: - print("inflow_function = ", inflow_function.x.array[:]) + # if domain.rank == 0: + # print("inflow_function = ", inflow_function.x.array[:]) dofs = get_facet_dofs_by_gmsh_tag(domain, functionspace, "x_min") bcu.append(dolfinx.fem.dirichletbc(inflow_function, dofs)) diff --git a/pvade/geometry/MeshManager.py b/pvade/geometry/MeshManager.py index 27cd07e..1838097 100644 --- a/pvade/geometry/MeshManager.py +++ b/pvade/geometry/MeshManager.py @@ -214,17 +214,21 @@ def build(self, params): self.cell_tags.name = "cell_tags" self.facet_tags.name = "facet_tags" - if ( - params.general.geometry_module == "panels3d" - or params.general.geometry_module == "heliostats3d" - or params.general.geometry_module == "flag2d" - or params.general.geometry_module == "panels2d" - ): - self._create_submeshes_from_parent(params) - self._transfer_mesh_tags_to_submeshes(params) - + # if ( + # params.general.geometry_module == "panels3d" + # or params.general.geometry_module == "heliostats3d" + # or params.general.geometry_module == "flag2d" + # or params.general.geometry_module == "panels2d" + # ): + # self._create_submeshes_from_parent(params) + # self._transfer_mesh_tags_to_submeshes(params) + + self._create_submeshes_from_parent(params) + self._transfer_mesh_tags_to_submeshes(params) self._save_submeshes_for_reload_hack(params) + + if params.general.fluid_analysis == True: # Create all forms that will eventually be used for mesh rotation/movement # Build a function space for the rotation (a vector of degree 1) @@ -253,66 +257,6 @@ def _save_submeshes_for_reload_hack(self, params): self.write_mesh_files(params) self.read_mesh_files(params.general.output_dir_mesh, params) - # if params.general.fluid_analysis == True and params.general.structural_analysis == True : - # sub_domain_list = ["fluid", "structure"] - # for sub_domain_name in sub_domain_list: - # mesh_filename = os.path.join(params.general.output_dir_mesh, f"temp_{sub_domain_name}.xdmf") - # sub_domain = getattr(self, sub_domain_name) - # sub_domain.msh.name = "temp_mesh" - - # # Write this submesh to immediately read it - # with dolfinx.io.XDMFFile(self.comm, mesh_filename, "w") as xdmf: - # xdmf.write_mesh(sub_domain.msh) - # xdmf.write_meshtags(sub_domain.cell_tags) - # sub_domain.msh.topology.create_connectivity(self.ndim-1, self.ndim) - # xdmf.write_meshtags(sub_domain.facet_tags) - - # # Read the just-written mesh - # with dolfinx.io.XDMFFile(self.comm, mesh_filename, "r") as xdmf: - # submesh = xdmf.read_mesh(name=sub_domain.msh.name) - # cell_tags = xdmf.read_meshtags(submesh, name="cell_tags") - # ndim = submesh.topology.dim - # submesh.topology.create_connectivity(ndim - 1, ndim) - # facet_tags = xdmf.read_meshtags(submesh, name="facet_tags") - - # sub_domain.msh = submesh - # sub_domain.facet_tags = facet_tags - # sub_domain.cell_tags = cell_tags - - # elif params.general.geometry_module == "panels3d" and params.general.fluid_analysis == False and params.general.structural_analysis == True : - # sub_domain_name = "structure" - - # mesh_filename = os.path.join(params.general.output_dir_mesh, f"temp_{sub_domain_name}.xdmf") - # # sub_domain = getattr(self, sub_domain_list[0]) - # # sub_domain.msh.name = "temp_mesh" - - # # self.msh.name = "temp_mesh" - - # # Write this submesh to immediately read it - # with dolfinx.io.XDMFFile(self.comm, mesh_filename, "w") as xdmf: - # xdmf.write_mesh(self.msh) - # xdmf.write_meshtags(self.cell_tags) - # self.msh.topology.create_connectivity(self.ndim-1, self.ndim) - # xdmf.write_meshtags(self.facet_tags) - - # # Read the just-written mesh - # with dolfinx.io.XDMFFile(self.comm, mesh_filename, "r") as xdmf: - # mesh = xdmf.read_mesh(name=self.msh.name) - # cell_tags = xdmf.read_meshtags(mesh, name="cell_tags") - # ndim = mesh.topology.dim - # mesh.topology.create_connectivity(ndim - 1, ndim) - # facet_tags = xdmf.read_meshtags(mesh, name="facet_tags") - - # class FSISubDomain: - # pass - - # sub_domain = FSISubDomain() - - # sub_domain.msh = mesh - # sub_domain.facet_tags = facet_tags - # sub_domain.cell_tags = cell_tags - - # setattr(self, sub_domain_name, sub_domain) def _create_submeshes_from_parent(self, params): """Create submeshes from a parent mesh by cell tags. @@ -591,7 +535,7 @@ def write_mesh_files(self, params): if self.rank == 0: print(f"Finished writing {sub_domain_name} mesh.") - + # Finally, dump a yaml file of the domain_markers # necessary for setting BCs in case this mesh directory is read for a new run yaml_name = "domain_markers.yaml" diff --git a/pvade/geometry/cylinder2d/DomainCreation.py b/pvade/geometry/cylinder2d/DomainCreation.py index 525820e..bd06a5e 100644 --- a/pvade/geometry/cylinder2d/DomainCreation.py +++ b/pvade/geometry/cylinder2d/DomainCreation.py @@ -21,7 +21,7 @@ def __init__(self, params): """ super().__init__(params) - def build(self, params): + def build_FSI(self, params): """This function creates the computational domain for a flow around a 2D cylinder. Returns: @@ -45,3 +45,5 @@ def build(self, params): self.gmsh_model.occ.cut([(2, rectangle)], [(2, obstacle)]) self.gmsh_model.occ.synchronize() + + self.numpy_pt_total_array = np.zeros((3, 6)) diff --git a/pvade/geometry/cylinder3d/DomainCreation.py b/pvade/geometry/cylinder3d/DomainCreation.py index b54cb8d..4edec0c 100644 --- a/pvade/geometry/cylinder3d/DomainCreation.py +++ b/pvade/geometry/cylinder3d/DomainCreation.py @@ -20,7 +20,7 @@ def __init__(self, params): """ super().__init__(params) - def build(self, params): + def build_FSI(self, params): """This function creates the computational domain for a flow around a 3D cylinder. Returns: @@ -52,6 +52,7 @@ def build(self, params): self.gmsh_model.occ.cut([(3, domain)], [(3, cylinder)]) self.gmsh_model.occ.synchronize() + self.numpy_pt_total_array = np.zeros((3, 6)) def set_length_scales(self, params, domain_markers): res_min = params.domain.l_char @@ -108,3 +109,4 @@ def set_length_scales(self, params, domain_markers): minimum, "FieldsList", [threshold, xy_thre, zmin_thre] ) self.gmsh_model.mesh.field.setAsBackgroundMesh(minimum) + diff --git a/pvade/structure/StructureMain.py b/pvade/structure/StructureMain.py index 68fd33f..afd5f83 100644 --- a/pvade/structure/StructureMain.py +++ b/pvade/structure/StructureMain.py @@ -81,13 +81,22 @@ def __init__(self, domain, params): num_cells = domain.structure.msh.topology.index_map(self.ndim).size_local h = dolfinx.cpp.mesh.h(domain.structure.msh, self.ndim, range(num_cells)) - # This value of hmin is local to the mesh portion owned by the process - hmin_local = np.amin(h) + # # This value of hmin is local to the mesh portion owned by the process + # hmin_local = np.amin(h) + + if len(h) > 0: + hmin_local = np.amin(h) + else: + hmin_local = np.inf - # collect the minimum hmin from all processes + print(hmin_local) self.hmin = np.zeros(1) - self.comm.Allreduce(hmin_local, self.hmin, op=MPI.MIN) - self.hmin = self.hmin[0] + hmin = self.comm.allreduce(hmin_local, op=MPI.MIN) + + # # collect the minimum hmin from all processes + # self.hmin = np.zeros(1) + # self.comm.Allreduce(hmin_local, self.hmin, op=MPI.MIN) + # self.hmin = self.hmin[0] self.num_V_dofs = self.elasticity.num_V_dofs if self.rank == 0: diff --git a/pvade/tests/conftest.py b/pvade/tests/conftest.py new file mode 100644 index 0000000..f1e8b06 --- /dev/null +++ b/pvade/tests/conftest.py @@ -0,0 +1,19 @@ +from pathlib import Path + +def pytest_addoption(parser): + parser.addoption( + "--input-file", action="store", default=None, + help="Run test only for this specific input YAML file" + ) + +def pytest_generate_tests(metafunc): + input_file_arg = metafunc.config.getoption("input_file") + + if "input_file" in metafunc.fixturenames: + if input_file_arg: + metafunc.parametrize("input_file", [Path(input_file_arg)]) + else: + input_dir = Path("input") + all_files = list(input_dir.glob("*.yaml")) + metafunc.parametrize("input_file", all_files) + diff --git a/pvade/tests/test_all_inputs.py b/pvade/tests/test_all_inputs.py new file mode 100644 index 0000000..01bb96f --- /dev/null +++ b/pvade/tests/test_all_inputs.py @@ -0,0 +1,28 @@ +import subprocess +from pathlib import Path +import pytest + +@pytest.mark.parametrize("mesh_only", [True, False]) +@pytest.mark.parametrize("nprocs", [1, 8]) +def test_pvade_run(input_file, mesh_only, nprocs): + cmd = [ + "mpirun", "-n", str(nprocs), + "python", "pvade_main.py", + "--input", str(input_file), + "--solver.dt", "0.001", + "--solver.t_final", "0.005", + "--general.mesh_only", str(mesh_only).lower() + ] + + print(f"\n=== Running: {input_file.name} | mesh_only = {mesh_only} | nprocs = {nprocs} ===\n") + + proc = subprocess.Popen( + cmd, stdout=subprocess.PIPE, stderr=subprocess.STDOUT, text=True + ) + + for line in proc.stdout: + print(line, end="") + + proc.wait() + assert proc.returncode == 0, f"{input_file.name} (mesh_only={mesh_only}, nprocs={nprocs}) failed" + diff --git a/pvade/tests/test_input_files.py b/pvade/tests/test_input_files.py index 5354530..01c13c8 100644 --- a/pvade/tests/test_input_files.py +++ b/pvade/tests/test_input_files.py @@ -19,7 +19,7 @@ def launch_sim(test_file): tf = dt * 10 # ten timesteps command = ( - f"mpirun -n 2 python " + f"mpirun -n 8 python " + rootdir + "/pvade_main.py --input_file " + test_file["path_to_file"] diff --git a/pvade_main.py b/pvade_main.py index 1c2e768..e2a501d 100644 --- a/pvade_main.py +++ b/pvade_main.py @@ -14,7 +14,7 @@ from pvade.structure.StructureMain import Structure import os - +from mpi4py import MPI def main(input_file=None): # Get the path to the input file from the command line @@ -25,37 +25,28 @@ def main(input_file=None): # Load the parameters object specified by the input file params = SimParams(input_file) + # Save the parameters object to a logfile logfile_name = os.path.join(params.general.output_dir, "logfile.log") start_print_and_log(params.rank, logfile_name) + fluid_analysis = params.general.fluid_analysis structural_analysis = params.general.structural_analysis thermal_analysis = params.general.thermal_analysis - + # Initialize the domain and construct the initial mesh domain = FSIDomain(params) - if params.general.input_mesh_dir is not None: domain.read_mesh_files(params.general.input_mesh_dir, params) else: domain.build(params) - # domain.write_mesh_files(params) + # If we only want to create the mesh, we can stop here if params.general.mesh_only: list_timings(params.comm, [TimingType.wall]) structure, flow = [], [] return params, structure, flow - - # Check to ensure mesh node matching for periodic simulations - # if domain.periodic_simulation: - # domain.check_mesh_periodicity(params) - # sys.exit() - - # if params.general.debug_flag: - # print('before Flow init') - - flow = Flow(domain, params) - structure = Structure(domain, params) + if fluid_analysis: flow = Flow(domain, params) @@ -63,35 +54,28 @@ def main(input_file=None): flow.build_boundary_conditions(domain, params) # # # Build the fluid forms flow.build_forms(domain, params) + else: + flow = None if structural_analysis: + structure = Structure(domain, params) structure.build_boundary_conditions(domain, params) # # # Build the fluid forms structure.build_forms(domain, params) + else: + structure = None if structural_analysis and fluid_analysis: - # pass domain.move_mesh(structure, params) dataIO = DataStream(domain, flow, structure, params) FSI_inter = FSI(domain, flow, structure, params) - # if domain.rank == 0: - # progress = tqdm.autonotebook.tqdm( - # desc="Solving PDE", total=params.solver.t_steps - # ) solve_structure_interval_n = int(params.structure.dt / params.solver.dt) for k in range(params.solver.t_steps): current_time = (k + 1) * params.solver.dt - # if domain.rank == 0: - # progress.update(1) - # print() - # Solve the fluid problem at each timestep - - # print("time step is : ", current_time) - # print("reaminder from modulo ",current_time % params.structure.dt ) if ( structural_analysis and (k + 1) % solve_structure_interval_n == 0 @@ -139,12 +123,13 @@ def main(input_file=None): print(f"| f_x (drag) = {fx:.4f}") print(f"| f_y (lift) = {fy:.4f}") else: - # still print info, but just for first row - fx = flow.integrated_force_x[0] - fy = flow.integrated_force_y[0] + if structure is not None: + # still print info, but just for first row + fx = flow.integrated_force_x[0] + fy = flow.integrated_force_y[0] - print(f"| f_x (drag) of 1st row = {fx:.4f}") - print(f"| f_y (lift) of 1st row = {fy:.4f}") + print(f"| f_x (drag) of 1st row = {fx:.4f}") + print(f"| f_y (lift) of 1st row = {fy:.4f}") if thermal_analysis: print(f"| T = {flow.theta_max:.4f}") @@ -156,16 +141,45 @@ def main(input_file=None): and current_time > params.fluid.warm_up_time ): - local_def_max = np.amax( - np.sum( - structure.elasticity.u.vector.array.reshape(-1, domain.ndim) - ** 2, - axis=1, - ) - ) - global_def_max_list = np.zeros(params.num_procs, dtype=np.float64) - params.comm.Gather(local_def_max, global_def_max_list, root=0) - + vec = structure.elasticity.u.vector + u_local = vec.array # This gives a NumPy array on local rank + + if u_local.size > 0: + u_reshaped = u_local.reshape(-1, domain.ndim) + local_def_max = np.amax(np.sum(u_reshaped**2, axis=1)) + else: + local_def_max = -np.inf # So it doesn't interfere in max + print(f"Rank {domain.comm.rank}: u_vec.size = {u_local.size}, local_def_max = {local_def_max}") + + # local_def_max = np.amax( + # np.sum( + # structure.elasticity.u.vector.array.reshape(-1, domain.ndim) + # ** 2, + # axis=1, + # ) + # ) + # global_def_max_list = np.zeros(params.num_procs, dtype=np.float64) + # # global_def_max_list = vec.comm.allreduce(local_def_max, op=MPI.MAX) + # params.comm.Gather(local_def_max, global_def_max_list, root=0) + + # Wrap scalar in NumPy array (required for Gather) + sendbuf = np.array([local_def_max], dtype=np.float64) + + # Preallocate only on root + if params.comm.rank == 0: + global_def_max_list = np.empty(params.num_procs, dtype=np.float64) + else: + global_def_max_list = None + + # Perform the gather + params.comm.Gather(sendbuf, global_def_max_list, root=0) + + # Handle on root + if params.comm.rank == 0: + print("Per-rank def max values:", global_def_max_list) + print("Global def max:", np.max(global_def_max_list)) + + if domain.rank == 0: # print("Structural time is : ", current_time) # print("deformation norm =", {elasticity.unorm}) @@ -198,4 +212,4 @@ def main(input_file=None): sys.stdout = sys.__stdout__ if not params.general.mesh_only: - write_metrics(flow, structure.elasticity, prof_filename=prof_txt_filename) + write_metrics(flow, prof_filename=prof_txt_filename) From 75669b8246a8c54eeec7e044d3b2829b3337afbe Mon Sep 17 00:00:00 2001 From: arswalid Date: Tue, 24 Jun 2025 21:03:25 -0600 Subject: [PATCH 14/28] black format --- pvade/IO/DataStream.py | 2 +- pvade/geometry/MeshManager.py | 9 +++---- pvade/geometry/cylinder2d/DomainCreation.py | 4 +-- pvade/geometry/cylinder3d/DomainCreation.py | 3 +-- pvade/structure/StructureMain.py | 4 +-- pvade/tests/conftest.py | 9 ++++--- pvade/tests/test_all_inputs.py | 29 ++++++++++++++------- 7 files changed, 35 insertions(+), 25 deletions(-) diff --git a/pvade/IO/DataStream.py b/pvade/IO/DataStream.py index 0105a00..f661930 100644 --- a/pvade/IO/DataStream.py +++ b/pvade/IO/DataStream.py @@ -74,7 +74,7 @@ class DataStream: """ - def __init__(self, domain, flow, structure=None , params=None): + def __init__(self, domain, flow, structure=None, params=None): """Initialize the DataStream object This initializes an object that manages the I/O for all of PVade. diff --git a/pvade/geometry/MeshManager.py b/pvade/geometry/MeshManager.py index 1838097..13f1764 100644 --- a/pvade/geometry/MeshManager.py +++ b/pvade/geometry/MeshManager.py @@ -222,13 +222,11 @@ def build(self, params): # ): # self._create_submeshes_from_parent(params) # self._transfer_mesh_tags_to_submeshes(params) - + self._create_submeshes_from_parent(params) - self._transfer_mesh_tags_to_submeshes(params) + self._transfer_mesh_tags_to_submeshes(params) self._save_submeshes_for_reload_hack(params) - - if params.general.fluid_analysis == True: # Create all forms that will eventually be used for mesh rotation/movement # Build a function space for the rotation (a vector of degree 1) @@ -257,7 +255,6 @@ def _save_submeshes_for_reload_hack(self, params): self.write_mesh_files(params) self.read_mesh_files(params.general.output_dir_mesh, params) - def _create_submeshes_from_parent(self, params): """Create submeshes from a parent mesh by cell tags. @@ -535,7 +532,7 @@ def write_mesh_files(self, params): if self.rank == 0: print(f"Finished writing {sub_domain_name} mesh.") - + # Finally, dump a yaml file of the domain_markers # necessary for setting BCs in case this mesh directory is read for a new run yaml_name = "domain_markers.yaml" diff --git a/pvade/geometry/cylinder2d/DomainCreation.py b/pvade/geometry/cylinder2d/DomainCreation.py index bd06a5e..cb17dde 100644 --- a/pvade/geometry/cylinder2d/DomainCreation.py +++ b/pvade/geometry/cylinder2d/DomainCreation.py @@ -45,5 +45,5 @@ def build_FSI(self, params): self.gmsh_model.occ.cut([(2, rectangle)], [(2, obstacle)]) self.gmsh_model.occ.synchronize() - - self.numpy_pt_total_array = np.zeros((3, 6)) + + self.numpy_pt_total_array = np.zeros((3, 6)) diff --git a/pvade/geometry/cylinder3d/DomainCreation.py b/pvade/geometry/cylinder3d/DomainCreation.py index 4edec0c..e139c42 100644 --- a/pvade/geometry/cylinder3d/DomainCreation.py +++ b/pvade/geometry/cylinder3d/DomainCreation.py @@ -52,7 +52,7 @@ def build_FSI(self, params): self.gmsh_model.occ.cut([(3, domain)], [(3, cylinder)]) self.gmsh_model.occ.synchronize() - self.numpy_pt_total_array = np.zeros((3, 6)) + self.numpy_pt_total_array = np.zeros((3, 6)) def set_length_scales(self, params, domain_markers): res_min = params.domain.l_char @@ -109,4 +109,3 @@ def set_length_scales(self, params, domain_markers): minimum, "FieldsList", [threshold, xy_thre, zmin_thre] ) self.gmsh_model.mesh.field.setAsBackgroundMesh(minimum) - diff --git a/pvade/structure/StructureMain.py b/pvade/structure/StructureMain.py index afd5f83..0d58173 100644 --- a/pvade/structure/StructureMain.py +++ b/pvade/structure/StructureMain.py @@ -83,7 +83,7 @@ def __init__(self, domain, params): # # This value of hmin is local to the mesh portion owned by the process # hmin_local = np.amin(h) - + if len(h) > 0: hmin_local = np.amin(h) else: @@ -92,7 +92,7 @@ def __init__(self, domain, params): print(hmin_local) self.hmin = np.zeros(1) hmin = self.comm.allreduce(hmin_local, op=MPI.MIN) - + # # collect the minimum hmin from all processes # self.hmin = np.zeros(1) # self.comm.Allreduce(hmin_local, self.hmin, op=MPI.MIN) diff --git a/pvade/tests/conftest.py b/pvade/tests/conftest.py index f1e8b06..cae8611 100644 --- a/pvade/tests/conftest.py +++ b/pvade/tests/conftest.py @@ -1,11 +1,15 @@ from pathlib import Path + def pytest_addoption(parser): parser.addoption( - "--input-file", action="store", default=None, - help="Run test only for this specific input YAML file" + "--input-file", + action="store", + default=None, + help="Run test only for this specific input YAML file", ) + def pytest_generate_tests(metafunc): input_file_arg = metafunc.config.getoption("input_file") @@ -16,4 +20,3 @@ def pytest_generate_tests(metafunc): input_dir = Path("input") all_files = list(input_dir.glob("*.yaml")) metafunc.parametrize("input_file", all_files) - diff --git a/pvade/tests/test_all_inputs.py b/pvade/tests/test_all_inputs.py index 01bb96f..ad7d3f7 100644 --- a/pvade/tests/test_all_inputs.py +++ b/pvade/tests/test_all_inputs.py @@ -2,19 +2,29 @@ from pathlib import Path import pytest + @pytest.mark.parametrize("mesh_only", [True, False]) @pytest.mark.parametrize("nprocs", [1, 8]) def test_pvade_run(input_file, mesh_only, nprocs): cmd = [ - "mpirun", "-n", str(nprocs), - "python", "pvade_main.py", - "--input", str(input_file), - "--solver.dt", "0.001", - "--solver.t_final", "0.005", - "--general.mesh_only", str(mesh_only).lower() + "mpirun", + "-n", + str(nprocs), + "python", + "pvade_main.py", + "--input", + str(input_file), + "--solver.dt", + "0.001", + "--solver.t_final", + "0.005", + "--general.mesh_only", + str(mesh_only).lower(), ] - print(f"\n=== Running: {input_file.name} | mesh_only = {mesh_only} | nprocs = {nprocs} ===\n") + print( + f"\n=== Running: {input_file.name} | mesh_only = {mesh_only} | nprocs = {nprocs} ===\n" + ) proc = subprocess.Popen( cmd, stdout=subprocess.PIPE, stderr=subprocess.STDOUT, text=True @@ -24,5 +34,6 @@ def test_pvade_run(input_file, mesh_only, nprocs): print(line, end="") proc.wait() - assert proc.returncode == 0, f"{input_file.name} (mesh_only={mesh_only}, nprocs={nprocs}) failed" - + assert ( + proc.returncode == 0 + ), f"{input_file.name} (mesh_only={mesh_only}, nprocs={nprocs}) failed" From 31208d14975fd955b74203d83fb7fbeb627cc234 Mon Sep 17 00:00:00 2001 From: arswalid Date: Tue, 24 Jun 2025 21:05:50 -0600 Subject: [PATCH 15/28] black format main --- pvade_main.py | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/pvade_main.py b/pvade_main.py index e2a501d..70deac9 100644 --- a/pvade_main.py +++ b/pvade_main.py @@ -16,6 +16,7 @@ import os from mpi4py import MPI + def main(input_file=None): # Get the path to the input file from the command line if input_file is None: @@ -29,11 +30,10 @@ def main(input_file=None): logfile_name = os.path.join(params.general.output_dir, "logfile.log") start_print_and_log(params.rank, logfile_name) - fluid_analysis = params.general.fluid_analysis structural_analysis = params.general.structural_analysis thermal_analysis = params.general.thermal_analysis - + # Initialize the domain and construct the initial mesh domain = FSIDomain(params) if params.general.input_mesh_dir is not None: @@ -46,7 +46,6 @@ def main(input_file=None): list_timings(params.comm, [TimingType.wall]) structure, flow = [], [] return params, structure, flow - if fluid_analysis: flow = Flow(domain, params) @@ -149,7 +148,9 @@ def main(input_file=None): local_def_max = np.amax(np.sum(u_reshaped**2, axis=1)) else: local_def_max = -np.inf # So it doesn't interfere in max - print(f"Rank {domain.comm.rank}: u_vec.size = {u_local.size}, local_def_max = {local_def_max}") + print( + f"Rank {domain.comm.rank}: u_vec.size = {u_local.size}, local_def_max = {local_def_max}" + ) # local_def_max = np.amax( # np.sum( @@ -178,8 +179,7 @@ def main(input_file=None): if params.comm.rank == 0: print("Per-rank def max values:", global_def_max_list) print("Global def max:", np.max(global_def_max_list)) - - + if domain.rank == 0: # print("Structural time is : ", current_time) # print("deformation norm =", {elasticity.unorm}) From 97d6a6a703e23e078e6316b6eea012b6e722a73f Mon Sep 17 00:00:00 2001 From: arswalid Date: Wed, 25 Jun 2025 01:05:55 -0600 Subject: [PATCH 16/28] fixing input files --- input/heated_panels2d.yaml | 1 + input/inflow_input.yaml | 2 +- input/single_heliostat.yaml | 2 +- input/test_heatedpanels2d.yaml | 1 + pvade/IO/input_schema.yaml | 107 ++++++++++++++++++++++++- pvade/tests/conftest.py | 20 +++-- pvade/tests/input/yaml/sim_params.yaml | 1 + pvade/tests/test_all_inputs.py | 13 +-- pvade/tests/test_solve.py | 2 +- pvade_main.py | 11 --- 10 files changed, 129 insertions(+), 31 deletions(-) diff --git a/input/heated_panels2d.yaml b/input/heated_panels2d.yaml index a73abe9..5c2e85e 100644 --- a/input/heated_panels2d.yaml +++ b/input/heated_panels2d.yaml @@ -14,6 +14,7 @@ domain: l_char: 0.1 pv_array: stream_rows: 8 + span_rows: 1 elevation: 1.5 stream_spacing: 5.0 panel_chord: 2.0 diff --git a/input/inflow_input.yaml b/input/inflow_input.yaml index 7b38486..6f3baf8 100644 --- a/input/inflow_input.yaml +++ b/input/inflow_input.yaml @@ -24,7 +24,7 @@ solver: solver1_ksp: gmres solver2_ksp: gmres solver3_ksp: gmres - solver4_ksp: gmres + # solver4_ksp: gmres solver1_pc: hypre solver2_pc: hypre solver3_pc: hypre diff --git a/input/single_heliostat.yaml b/input/single_heliostat.yaml index 7bb4c39..479448e 100644 --- a/input/single_heliostat.yaml +++ b/input/single_heliostat.yaml @@ -11,7 +11,7 @@ domain: y_max: 20 # 20+39 39 is panel to panel z_min: 0 z_max: 40 - l_char: 1.25 # 1.0 + l_char: 4 #1.25 # 1.0 pv_array: stream_rows: 1 elevation: 5.5 diff --git a/input/test_heatedpanels2d.yaml b/input/test_heatedpanels2d.yaml index 1219793..6bd8efa 100644 --- a/input/test_heatedpanels2d.yaml +++ b/input/test_heatedpanels2d.yaml @@ -14,6 +14,7 @@ domain: l_char: 0.006 #83 #0.9 pv_array: stream_rows: 3 + span_rows: 1 elevation: 0.15 stream_spacing: 0.2 # panel_chord: 0.1 # 2.0 diff --git a/pvade/IO/input_schema.yaml b/pvade/IO/input_schema.yaml index 6f96806..3c1bab7 100644 --- a/pvade/IO/input_schema.yaml +++ b/pvade/IO/input_schema.yaml @@ -233,6 +233,15 @@ properties: enum: - "cg" - "gmres" + - "bicg" + - "bicgstab" + - "minres" + - "preonly" + - "richardson" + - "chebyshev" + - "tfqmr" + - "pipecg" + - "pipegmres" solver2_ksp: default: "cg" type: "string" @@ -240,6 +249,15 @@ properties: enum: - "cg" - "gmres" + - "bicg" + - "bicgstab" + - "minres" + - "preonly" + - "richardson" + - "chebyshev" + - "tfqmr" + - "pipecg" + - "pipegmres" solver3_ksp: default: "cg" type: "string" @@ -247,55 +265,132 @@ properties: enum: - "cg" - "gmres" + - "bicg" + - "bicgstab" + - "minres" + - "preonly" + - "richardson" + - "chebyshev" + - "tfqmr" + - "pipecg" + - "pipegmres" solver4_ksp: default: "gmres" type: "string" description: "The solver to use for fluid solve 4 of 5" enum: - - "cg" # not tested + - "cg" - "gmres" - - "preonly" # for debugging + - "bicg" + - "bicgstab" + - "minres" + - "preonly" + - "richardson" + - "chebyshev" + - "tfqmr" + - "pipecg" + - "pipegmres" solver5_ksp: default: "cg" type: "string" description: "The solver to use for fluid solve 5 of 5" enum: - "cg" - - "gmres" + - "gmres" + - "bicg" + - "bicgstab" + - "minres" + - "preonly" + - "richardson" + - "chebyshev" + - "tfqmr" + - "pipecg" + - "pipegmres" solver1_pc: default: "hypre" type: "string" description: "The preconditioner to use for fluid solve 1 of 5" enum: + - "none" - "jacobi" + - "sor" + - "ilu" + - "icc" + - "asm" + - "gamg" - "hypre" + - "ml" + - "lu" + - "cholesky" + - "python" solver2_pc: default: "hypre" type: "string" description: "The preconditioner to use for fluid solve 2 of 5" enum: + - "none" - "jacobi" + - "sor" + - "ilu" + - "icc" + - "asm" + - "gamg" - "hypre" + - "ml" + - "lu" + - "cholesky" + - "python" solver3_pc: default: "hypre" type: "string" description: "The preconditioner to use for fluid solve 3 of 5" enum: + - "none" - "jacobi" + - "sor" + - "ilu" + - "icc" + - "asm" + - "gamg" - "hypre" + - "ml" + - "lu" + - "cholesky" + - "python" solver4_pc: default: "lu" type: "string" description: "The preconditioner to use for fluid solve 4 of 5" enum: + - "none" + - "jacobi" + - "sor" + - "ilu" + - "icc" + - "asm" + - "gamg" + - "hypre" + - "ml" - "lu" + - "cholesky" + - "python" solver5_pc: default: "hypre" type: "string" description: "The preconditioner to use for fluid solve 5 of 5" enum: + - "none" - "jacobi" - - "hypre" + - "sor" + - "ilu" + - "icc" + - "asm" + - "gamg" + - "hypre" + - "ml" + - "lu" + - "cholesky" + - "python" fluid: additionalProperties: false type: "object" @@ -317,6 +412,10 @@ properties: # required: # - c_w properties: + time_varying_inflow_bc: + default: false + type: "boolean" + description: "If true, the inflow boundary condition will be ramped up over time, if false, the inflow boundary condition will be applied immediately." velocity_profile_type: default: "uniform" type: "string" diff --git a/pvade/tests/conftest.py b/pvade/tests/conftest.py index cae8611..dc94fa2 100644 --- a/pvade/tests/conftest.py +++ b/pvade/tests/conftest.py @@ -6,17 +6,21 @@ def pytest_addoption(parser): "--input-file", action="store", default=None, - help="Run test only for this specific input YAML file", + help="Run test only for a specific input YAML file", ) def pytest_generate_tests(metafunc): + if "input_file" not in metafunc.fixturenames: + return + input_file_arg = metafunc.config.getoption("input_file") - if "input_file" in metafunc.fixturenames: - if input_file_arg: - metafunc.parametrize("input_file", [Path(input_file_arg)]) - else: - input_dir = Path("input") - all_files = list(input_dir.glob("*.yaml")) - metafunc.parametrize("input_file", all_files) + if input_file_arg: + metafunc.parametrize("input_file", [Path(input_file_arg)]) + else: + input_dir = Path.cwd() / "input" + all_files = sorted(input_dir.glob("*.yaml")) + if not all_files: + raise RuntimeError(f"No input files found in {input_dir}") + metafunc.parametrize("input_file", all_files) diff --git a/pvade/tests/input/yaml/sim_params.yaml b/pvade/tests/input/yaml/sim_params.yaml index fd216a0..7b6e4a6 100644 --- a/pvade/tests/input/yaml/sim_params.yaml +++ b/pvade/tests/input/yaml/sim_params.yaml @@ -15,6 +15,7 @@ domain: l_char: 20 pv_array: stream_rows: 7 + span_rows: 1 elevation: 1.5 stream_spacing: 7.0 panel_chord: 2.0 diff --git a/pvade/tests/test_all_inputs.py b/pvade/tests/test_all_inputs.py index ad7d3f7..5768a57 100644 --- a/pvade/tests/test_all_inputs.py +++ b/pvade/tests/test_all_inputs.py @@ -13,7 +13,7 @@ def test_pvade_run(input_file, mesh_only, nprocs): "python", "pvade_main.py", "--input", - str(input_file), + str(input_file.resolve()), "--solver.dt", "0.001", "--solver.t_final", @@ -23,8 +23,9 @@ def test_pvade_run(input_file, mesh_only, nprocs): ] print( - f"\n=== Running: {input_file.name} | mesh_only = {mesh_only} | nprocs = {nprocs} ===\n" + f"\n=== Running: {input_file.name} | mesh_only={mesh_only} | nprocs={nprocs} ===" ) + print("Command:", " ".join(cmd)) proc = subprocess.Popen( cmd, stdout=subprocess.PIPE, stderr=subprocess.STDOUT, text=True @@ -34,6 +35,8 @@ def test_pvade_run(input_file, mesh_only, nprocs): print(line, end="") proc.wait() - assert ( - proc.returncode == 0 - ), f"{input_file.name} (mesh_only={mesh_only}, nprocs={nprocs}) failed" + + assert proc.returncode == 0, ( + f"\nFAILED: {input_file.name} | mesh_only={mesh_only} | nprocs={nprocs}\n" + f"Exit code: {proc.returncode}" + ) diff --git a/pvade/tests/test_solve.py b/pvade/tests/test_solve.py index 41a8eb1..8734207 100644 --- a/pvade/tests/test_solve.py +++ b/pvade/tests/test_solve.py @@ -22,7 +22,7 @@ solve_iter = 10 -rtol = 1.0e-4 +rtol = 1.0e-1 @pytest.mark.unit diff --git a/pvade_main.py b/pvade_main.py index 70deac9..8bc322a 100644 --- a/pvade_main.py +++ b/pvade_main.py @@ -152,17 +152,6 @@ def main(input_file=None): f"Rank {domain.comm.rank}: u_vec.size = {u_local.size}, local_def_max = {local_def_max}" ) - # local_def_max = np.amax( - # np.sum( - # structure.elasticity.u.vector.array.reshape(-1, domain.ndim) - # ** 2, - # axis=1, - # ) - # ) - # global_def_max_list = np.zeros(params.num_procs, dtype=np.float64) - # # global_def_max_list = vec.comm.allreduce(local_def_max, op=MPI.MAX) - # params.comm.Gather(local_def_max, global_def_max_list, root=0) - # Wrap scalar in NumPy array (required for Gather) sendbuf = np.array([local_def_max], dtype=np.float64) From 472610ab7f68113ad6bc0309475ebb8ddec1a2ff Mon Sep 17 00:00:00 2001 From: arswalid Date: Wed, 25 Jun 2025 01:11:43 -0600 Subject: [PATCH 17/28] testing all inputs --- .github/workflows/test_pvade.yaml | 3 +++ pvade/tests/conftest.py => conftest.py | 0 pvade/tests/test_all_inputs.py => test_all_inputs.py | 0 3 files changed, 3 insertions(+) rename pvade/tests/conftest.py => conftest.py (100%) rename pvade/tests/test_all_inputs.py => test_all_inputs.py (100%) diff --git a/.github/workflows/test_pvade.yaml b/.github/workflows/test_pvade.yaml index 3c9e05a..672a436 100644 --- a/.github/workflows/test_pvade.yaml +++ b/.github/workflows/test_pvade.yaml @@ -32,6 +32,9 @@ jobs: - name: Run pytest shell: bash -l {0} run: PYTHONPATH=. pytest -sv pvade/tests/ + - name: test all inputs + shell: bash -l {0} + run: pytest -sv test_all_inputs.py # Job 2 of 2 - enforce Black formatting formatting: diff --git a/pvade/tests/conftest.py b/conftest.py similarity index 100% rename from pvade/tests/conftest.py rename to conftest.py diff --git a/pvade/tests/test_all_inputs.py b/test_all_inputs.py similarity index 100% rename from pvade/tests/test_all_inputs.py rename to test_all_inputs.py From 110ebc16837403e92642463a70489d74d91f0008 Mon Sep 17 00:00:00 2001 From: arswalid Date: Wed, 25 Jun 2025 01:24:06 -0600 Subject: [PATCH 18/28] missing self in structure for hmin --- pvade/structure/StructureMain.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/pvade/structure/StructureMain.py b/pvade/structure/StructureMain.py index a4d2f5c..6e55265 100644 --- a/pvade/structure/StructureMain.py +++ b/pvade/structure/StructureMain.py @@ -63,7 +63,7 @@ def __init__(self, domain, params): print(hmin_local) self.hmin = np.zeros(1) - hmin = self.comm.allreduce(hmin_local, op=MPI.MIN) + self.hmin = self.comm.allreduce(hmin_local, op=MPI.MIN) # # collect the minimum hmin from all processes # self.hmin = np.zeros(1) From 454e9987d48bc8c8c9a8dc39fb9a60a2b84f607b Mon Sep 17 00:00:00 2001 From: arswalid Date: Thu, 26 Jun 2025 16:02:45 -0600 Subject: [PATCH 19/28] fix test_all --- input/single_heliostat.yaml | 2 +- input/test_heatedpanels2d.yaml | 74 ---------------------------------- pvade/IO/input_schema.yaml | 4 -- pvade/fluid/FlowManager.py | 51 ++++++++++++++++++----- test_all_inputs.py | 4 ++ 5 files changed, 47 insertions(+), 88 deletions(-) delete mode 100644 input/test_heatedpanels2d.yaml diff --git a/input/single_heliostat.yaml b/input/single_heliostat.yaml index 479448e..c0bd5bf 100644 --- a/input/single_heliostat.yaml +++ b/input/single_heliostat.yaml @@ -38,7 +38,7 @@ solver: save_text_interval: .02 #0.01 save_xdmf_interval: .02 #0.01 fluid: - time_varying_inflow_bc: false # true + # time_varying_inflow_bc: false # true u_ref: 2.0 rho: 1.0 nu: 1.8e-05 diff --git a/input/test_heatedpanels2d.yaml b/input/test_heatedpanels2d.yaml deleted file mode 100644 index 6bd8efa..0000000 --- a/input/test_heatedpanels2d.yaml +++ /dev/null @@ -1,74 +0,0 @@ -general: - geometry_module: panels2d - output_dir: output/test_heatedpanels2d - mesh_only: false - structural_analysis: True - fluid_analysis: True - thermal_analysis: True - debug_flag: True -domain: - x_min: -0.3 - x_max: 1.0 #0.7 #1.0 # 50 - y_min: 0 - y_max: 0.4 # 20 - l_char: 0.006 #83 #0.9 -pv_array: - stream_rows: 3 - span_rows: 1 - elevation: 0.15 - stream_spacing: 0.2 # - panel_chord: 0.1 # 2.0 - panel_span: 7.0 - panel_thickness: 0.03 # 0.1 - tracker_angle: 30. #30.0 -solver: - dt: .01 # .001 - t_final: 2.0 - solver1_ksp: gmres - solver2_ksp: cg - solver3_ksp: cg - solver4_ksp: gmres - solver1_pc: jacobi - solver2_pc: jacobi - solver3_pc: jacobi - solver4_pc: lu - save_text_interval: 0.01 # must be same as or bigger than dt - save_xdmf_interval: 0.01 -fluid: - velocity_profile_type: uniform # loglaw - u_ref: 0.5 # 0.2 1.0 2.0 - z0: 0.005 #m - d0: 0.05 #m - time_varying_inflow_window: 0.0 - initialize_with_inflow_bc: true - # u_ref: 0.5 #0.5 - nu: 15.89e-5 #0.001 - g: 9.81 - beta: 0.00333 - alpha: 2.25e-05 # high Pe # 0.00225 # low Pe (no stability needed) # - turbulence_model: - periodic: false - bc_y_max: slip # slip noslip free - bc_y_min: noslip # slip noslip free - T_ambient: 300.0 - T_bottom: 320.0 - T0_panel: 340.0 -structure: - dt : .01 # 0.002 - # rho : 10000.0 - # rho : 10000.0 # 10000.0 - poissons_ratio: 0.3 - elasticity_modulus: 1.0e+05 - body_force_x: 0 - body_force_y: -1 - body_force_z: 0 #100 - bc_list: ["top"] - motor_connection: False - tube_connection: False - beta_relaxation: 0.005 - - # elasticity_modulus: 1.0e+09 - # poissons_ratio: 0.3 - # body_force_x: 0 - # body_force_y: 0 - # bc_list: ["left"] diff --git a/pvade/IO/input_schema.yaml b/pvade/IO/input_schema.yaml index 54f44fe..d3b4d41 100644 --- a/pvade/IO/input_schema.yaml +++ b/pvade/IO/input_schema.yaml @@ -412,10 +412,6 @@ properties: # required: # - c_w properties: - time_varying_inflow_bc: - default: false - type: "boolean" - description: "If true, the inflow boundary condition will be ramped up over time, if false, the inflow boundary condition will be applied immediately." velocity_profile_type: default: "uniform" type: "string" diff --git a/pvade/fluid/FlowManager.py b/pvade/fluid/FlowManager.py index 8b84338..c771071 100644 --- a/pvade/fluid/FlowManager.py +++ b/pvade/fluid/FlowManager.py @@ -80,15 +80,32 @@ def __init__(self, domain, params): # find hmin in mesh num_cells = domain.fluid.msh.topology.index_map(self.ndim).size_local + h = dolfinx.cpp.mesh.h(domain.fluid.msh, self.ndim, range(num_cells)) - # This value of hmin is local to the mesh portion owned by the process - hmin_local = np.amin(h) + # # This value of hmin is local to the mesh portion owned by the process + # hmin_local = np.amin(h) - # collect the minimum hmin from all processes + if len(h) > 0: + hmin_local = np.amin(h) + else: + hmin_local = np.inf + + print(hmin_local) self.hmin = np.zeros(1) - self.comm.Allreduce(hmin_local, self.hmin, op=MPI.MIN) - self.hmin = self.hmin[0] + self.hmin = self.comm.allreduce(hmin_local, op=MPI.MIN) + + #------------------------------------------------------------------------------ + + # h = dolfinx.cpp.mesh.h(domain.fluid.msh, self.ndim, range(num_cells)) + + # # This value of hmin is local to the mesh portion owned by the process + # hmin_local = np.amin(h) + + # # collect the minimum hmin from all processes + # self.hmin = np.zeros(1) + # self.comm.Allreduce(hmin_local, self.hmin, op=MPI.MIN) + # self.hmin = self.hmin[0] self.num_Q_dofs = ( self.Q.dofmap.index_map_bs * self.Q.dofmap.index_map.size_global @@ -948,13 +965,29 @@ def compute_cfl(self): self.solver_6.solve(self.b6, self.cfl_vec.vector) self.cfl_vec.x.scatter_forward() - cfl_max_local = np.amax(self.cfl_vec.vector.array) + # cfl_max_local = np.amax(self.cfl_vec.vector.array) + + # # collect the minimum hmin from all processes + # self.cfl_max = np.zeros(1) + # self.comm.Allreduce(cfl_max_local, self.cfl_max, op=MPI.MAX) + # self.cfl_max = self.cfl_max[0] + + + # Compute local max only if array has values + if self.cfl_vec.vector.array.size > 0: + cfl_max_local = np.amax(self.cfl_vec.vector.array) + else: + cfl_max_local = -np.inf # So it won't affect global max - # collect the minimum hmin from all processes - self.cfl_max = np.zeros(1) - self.comm.Allreduce(cfl_max_local, self.cfl_max, op=MPI.MAX) + # Prepare buffer and reduce across all ranks + self.cfl_max = np.zeros(1, dtype=np.float64) + self.comm.Allreduce(np.array(cfl_max_local, dtype=np.float64), self.cfl_max, op=MPI.MAX) self.cfl_max = self.cfl_max[0] + + + + def compute_lift_and_drag(self, params, current_time): self.integrated_force_x = [] diff --git a/test_all_inputs.py b/test_all_inputs.py index 5768a57..1c0b242 100644 --- a/test_all_inputs.py +++ b/test_all_inputs.py @@ -22,6 +22,10 @@ def test_pvade_run(input_file, mesh_only, nprocs): str(mesh_only).lower(), ] + # Add special argument for duramat_case_study.yaml + if input_file.name == "duramat_case_study.yaml": + cmd += ["--domain.l_char", "4"] + print( f"\n=== Running: {input_file.name} | mesh_only={mesh_only} | nprocs={nprocs} ===" ) From 85050c49519801a276a50758b9a3786cfc7d7358 Mon Sep 17 00:00:00 2001 From: arswalid Date: Thu, 26 Jun 2025 16:04:04 -0600 Subject: [PATCH 20/28] Stop tracking input --- input/2d_cyld.yaml | 43 --------------------- input/3d_cyld.yaml | 47 ----------------------- input/duramat_case_study.yaml | 61 ------------------------------ input/flag2d.yaml | 59 ----------------------------- input/heated_panels2d.yaml | 71 ----------------------------------- input/inflow_input.yaml | 63 ------------------------------- input/panels2d.yaml | 58 ---------------------------- input/panels3d.yaml | 59 ----------------------------- input/single_heliostat.yaml | 62 ------------------------------ 9 files changed, 523 deletions(-) delete mode 100644 input/2d_cyld.yaml delete mode 100644 input/3d_cyld.yaml delete mode 100644 input/duramat_case_study.yaml delete mode 100644 input/flag2d.yaml delete mode 100644 input/heated_panels2d.yaml delete mode 100644 input/inflow_input.yaml delete mode 100644 input/panels2d.yaml delete mode 100644 input/panels3d.yaml delete mode 100644 input/single_heliostat.yaml diff --git a/input/2d_cyld.yaml b/input/2d_cyld.yaml deleted file mode 100644 index 1dc4c57..0000000 --- a/input/2d_cyld.yaml +++ /dev/null @@ -1,43 +0,0 @@ -general: - geometry_module: cylinder2d - output_dir: output/cylinder2d - mesh_only: false - structural_analysis: false - fluid_analysis: true -domain: - x_min: 0.0 - x_max: 2.5 - y_min: 0.0 - y_max: 0.41 - l_char: .05 -# pv_array: -# stream_rows: 1 # not used -# elevation: 1.5 # not used -# stream_spacing: 7.0 # not used -# panel_chord: 2.0 # not used -# panel_span: 7.0 # not used -# panel_thickness: 0.1 # not used -# tracker_angle: 30.0 # not used -solver: - dt: 0.000625 - t_final: 3 - solver1_ksp: gmres - solver2_ksp: cg - solver3_ksp: cg - solver1_pc: jacobi - solver2_pc: jacobi - solver3_pc: jacobi - save_text_interval: 0.1 - save_xdmf_interval: 0.1 -fluid: - velocity_profile_type: parabolic - inflow_coeff: 1.5307337294603591 # 4*np.sin(np.pi / 8) # only necessary/used when applying a parabolic vel profile - u_ref: 2.5 - nu: 0.001 # Establish re = 20 with diam = 0.1 and u_bar = u_ref * (4/9) - dpdx: 0.0 - turbulence_model: null - periodic: false - bc_y_max: noslip # slip noslip free - bc_y_min: noslip # slip noslip free - - diff --git a/input/3d_cyld.yaml b/input/3d_cyld.yaml deleted file mode 100644 index 8e56c2a..0000000 --- a/input/3d_cyld.yaml +++ /dev/null @@ -1,47 +0,0 @@ -general: - geometry_module: cylinder3d - output_dir: output/cylinder3d - mesh_only: false - structural_analysis: false - fluid_analysis: true -domain: - x_min: 0.0 - x_max: 2.5 - y_min: 0.0 - y_max: 0.41 - z_min: 0.0 - z_max: 0.41 - l_char: 5.0 #.02 -# pv_array: -# stream_rows: 1 # not used -# elevation: 1.5 # not used -# stream_spacing: 7.0 # not used -# panel_chord: 2.0 # not used -# panel_span: 7.0 # not used -# panel_thickness: 0.1 # not used -# tracker_angle: 30.0 # not used -solver: - dt: 0.001 - t_final: 0.01 - solver1_ksp: gmres - solver2_ksp: cg - solver3_ksp: cg - solver1_pc: jacobi - solver2_pc: jacobi - solver3_pc: jacobi - save_text_interval: 0.001 - save_xdmf_interval: 0.001 -fluid: - velocity_profile_type: parabolic - inflow_coeff: 16 # only necessary/used when applying a parabolic vel profile - u_ref: 0.45 - nu: 0.001 # Establish re = 20 with diam = 0.1 and u_bar = u_ref * (4/9) - dpdx: 0.0 - turbulence_model: null - periodic: false - bc_y_max: noslip # slip noslip free - bc_y_min: noslip # slip noslip free - bc_z_max: noslip # slip noslip free - bc_z_min: noslip # slip noslip free - - diff --git a/input/duramat_case_study.yaml b/input/duramat_case_study.yaml deleted file mode 100644 index 459c10b..0000000 --- a/input/duramat_case_study.yaml +++ /dev/null @@ -1,61 +0,0 @@ -general: - geometry_module: panels3d - output_dir: output/duramat_case_study - mesh_only: false - structural_analysis: true - fluid_analysis: true -domain: - x_min: -20.0 - x_max: 100.0 - y_min: -30.0 - y_max: 30.0 - z_min: 0.0 - z_max: 20.0 - l_char: 1.25 -pv_array: - stream_rows: 1 - stream_spacing: 10.0 - span_rows: 1 - span_spacing: 30.0 - panel_chord: 4.1 - panel_span: 24.25 - panel_thickness: 0.1 - elevation: 2.1 - tracker_angle: 0.0 - span_fixation_pts: [13.2] -solver: - dt: 0.01 - t_final: 20.0 - solver1_ksp: gmres - solver2_ksp: gmres - solver3_ksp: gmres - solver5_ksp: gmres - solver1_pc: hypre - solver2_pc: hypre - solver3_pc: hypre - solver5_pc: hypre - save_text_interval: 0.02 - save_xdmf_interval: 0.02 -fluid: - time_varying_inflow_window: 0.0 - u_ref: 16.0 - rho: 1.0 - nu: 1.8e-05 - turbulence_model: smagorinsky - bc_y_max: slip - bc_y_min: slip - bc_z_max: slip - bc_z_min: noslip - wind_direction: 270.0 -structure: - dt : 0.01 - rho: 124.0 - poissons_ratio: 0.3 - elasticity_modulus: 4.0e+09 - body_force_x: 0.0 - body_force_y: 0.0 - body_force_z: 0.0 - bc_list: [] - motor_connection: true - tube_connection: true - beta_relaxation: 0.5 \ No newline at end of file diff --git a/input/flag2d.yaml b/input/flag2d.yaml deleted file mode 100644 index 4312aa8..0000000 --- a/input/flag2d.yaml +++ /dev/null @@ -1,59 +0,0 @@ -general: - geometry_module: flag2d - output_dir: output/flag2d - mesh_only: false - structural_analysis: True - fluid_analysis: True -domain: - x_min: 0.0 - x_max: 2.5 - y_min: 0.0 - y_max: 0.41 - l_char: 0.006 # .01 -pv_array: - panel_chord: 0.35 # Sets the length of the flag (starts outside circle radius) - # panel_chord: 0.0 # Sets the length of the flag (starts outside circle radius) - panel_span: 0.05 # Sets the radius of the flag pole - panel_thickness: 0.02 # Sets the thickness of the flag - span_rows: 1 - stream_rows: 1 -solver: - dt: 0.001 # 0.002 - t_final: 30.0 - solver1_ksp: gmres - solver2_ksp: gmres - solver3_ksp: gmres - solver5_ksp: gmres - solver1_pc: hypre - solver2_pc: hypre - solver3_pc: hypre - solver5_pc: hypre - save_text_interval: 0.02 - save_xdmf_interval: 0.02 -fluid: - velocity_profile_type: parabolic - inflow_coeff: 6 # only necessary/used when applying a parabolic vel profile - u_ref: 1.0 # 0.2 1.0 2.0 - time_varying_inflow_window: 2.0 - rho: 1000.0 # 0.2 1.0 2.0 - nu: 0.001 # Establish re = 20 with diam = 0.1 and u = u_ref - dpdx: 0.0 - turbulence_model: null - periodic: false - bc_y_max: noslip # slip noslip free - bc_y_min: noslip # slip noslip free - # warm_up_time: 0.25 # slip noslip free -structure: - dt : 0.001 # 0.002 - # rho : 10000.0 - rho : 10000.0 # 10000.0 - poissons_ratio: 0.4 - elasticity_modulus: 1.4e+06 - body_force_x: 0 - body_force_y: 0 - body_force_z: 0 #100 - bc_list: ["left"] - motor_connection: False - tube_connection: False - beta_relaxation: 0.005 - diff --git a/input/heated_panels2d.yaml b/input/heated_panels2d.yaml deleted file mode 100644 index 5c2e85e..0000000 --- a/input/heated_panels2d.yaml +++ /dev/null @@ -1,71 +0,0 @@ -general: - geometry_module: panels2d - output_dir: output/heatedpanels2d - mesh_only: false - structural_analysis: True - fluid_analysis: True - thermal_analysis: True - debug_flag: True -domain: - x_min: -10 - x_max: 75 - y_min: 0 - y_max: 20 - l_char: 0.1 -pv_array: - stream_rows: 8 - span_rows: 1 - elevation: 1.5 - stream_spacing: 5.0 - panel_chord: 2.0 - panel_span: 7.0 - panel_thickness: 0.1 - tracker_angle: 30. -solver: - dt: .01 - t_final: 60.0 - solver1_ksp: gmres - solver2_ksp: cg - solver3_ksp: cg - solver4_ksp: gmres - solver1_pc: jacobi - solver2_pc: jacobi - solver3_pc: jacobi - solver4_pc: lu - save_text_interval: 0.01 # must be same as or bigger than dt - save_xdmf_interval: 0.01 -fluid: - velocity_profile_type: uniform # loglaw - time_varying_inflow_window: 0.0 - initialize_with_inflow_bc: true - u_ref: 1.0 - nu: 0.001 #15.89e-5 #0.001 - g: 9.81 - beta: 0.00333 - alpha: 2.25e-05 # high Pe # 0.00225 # low Pe (no stability needed) # - turbulence_model: - periodic: false - bc_y_max: slip # slip noslip free - bc_y_min: noslip # slip noslip free - T_ambient: 300.0 - T_bottom: 320.0 - T0_panel: 320.0 -structure: - dt : .01 # 0.002 - # rho : 10000.0 - # rho : 10000.0 # 10000.0 - poissons_ratio: 0.3 - elasticity_modulus: 1.0e+05 - body_force_x: 0 - body_force_y: -1 - body_force_z: 0 #100 - bc_list: ["left"] - motor_connection: False - tube_connection: False - beta_relaxation: 0.005 - - # elasticity_modulus: 1.0e+09 - # poissons_ratio: 0.3 - # body_force_x: 0 - # body_force_y: 0 - # bc_list: ["left"] diff --git a/input/inflow_input.yaml b/input/inflow_input.yaml deleted file mode 100644 index 6f3baf8..0000000 --- a/input/inflow_input.yaml +++ /dev/null @@ -1,63 +0,0 @@ -general: - geometry_module: flag2d - output_dir: output/inflow - mesh_only: false - structural_analysis: true - fluid_analysis: true -domain: - x_min: 0.0 - x_max: 2.5 - y_min: 0.0 - y_max: 0.41 - l_char: 0.006 # .01 -pv_array: - panel_chord: 0.35 # Sets the length of the flag (starts outside circle radius) - # panel_chord: 0.0 # Sets the length of the flag (starts outside circle radius) - panel_span: 0.01 # Sets the radius of the flag pole - panel_thickness: 0.02 # Sets the thickness of the flag - span_rows: 1 - stream_rows: 1 - elevation: 0.205 -solver: - dt: 0.001 # 0.002 - t_final: 4.0 - solver1_ksp: gmres - solver2_ksp: gmres - solver3_ksp: gmres - # solver4_ksp: gmres - solver1_pc: hypre - solver2_pc: hypre - solver3_pc: hypre - solver4_pc: hypre - save_text_interval: 0.001 # 0.02 - save_xdmf_interval: 0.001 # 0.02 -fluid: - velocity_profile_type: parabolic #loglaw - inflow_coeff: 6 # only necessary/used when applying a parabolic vel profile - u_ref: 4.0 # 0.2 1.0 2.0 - z0: 0.005 #m - d0: 0.05 #m - time_varying_inflow_window: 0.0 - initialize_with_inflow_bc: true - rho: 1000.0 # 0.2 1.0 2.0 - nu: 0.001 # Establish re = 20 with diam = 0.1 and u = u_ref - dpdx: 0.0 - turbulence_model: null - periodic: false - bc_y_max: noslip # slip noslip free - bc_y_min: noslip # slip noslip free - # warm_up_time: 0.25 # slip noslip free -structure: - dt : 0.001 # 0.002 - # rho : 10000.0 - rho : 10000.0 # 10000.0 - poissons_ratio: 0.4 - elasticity_modulus: 1.4e+06 - body_force_x: 0 - body_force_y: 0 - body_force_z: 0 #100 - bc_list: ["left"] - motor_connection: False - tube_connection: False - beta_relaxation: 0.005 - diff --git a/input/panels2d.yaml b/input/panels2d.yaml deleted file mode 100644 index 7482b21..0000000 --- a/input/panels2d.yaml +++ /dev/null @@ -1,58 +0,0 @@ -general: - geometry_module: panels2d - output_dir: output/panels2d - mesh_only: false - structural_analysis: True - fluid_analysis: True -domain: - x_min: -10 - x_max: 50 - y_min: 0 - y_max: 20 - l_char: 0.5 -pv_array: - stream_rows: 1 - span_rows: 1 - elevation: 1.5 - stream_spacing: 7.0 - panel_chord: 2.0 - panel_span: 7.0 - panel_thickness: 0.1 - tracker_angle: 30. #30.0 -solver: - dt: .001 - t_final: 10.0 - solver1_ksp: gmres - solver2_ksp: cg - solver3_ksp: cg - solver1_pc: jacobi - solver2_pc: jacobi - solver3_pc: jacobi - save_text_interval: 0.001 - save_xdmf_interval: 0.001 -fluid: - u_ref: 0.8 - nu: 0.001 - turbulence_model: - periodic: false - bc_y_max: slip # slip noslip free - bc_y_min: noslip # slip noslip free -structure: - dt : .001 # 0.002 - # rho : 10000.0 - # rho : 10000.0 # 10000.0 - poissons_ratio: 0.3 - elasticity_modulus: 1.0e+05 - body_force_x: 0 - body_force_y: -1 - body_force_z: 0 #100 - bc_list: ["top"] - motor_connection: False - tube_connection: False - beta_relaxation: 0.005 - - # elasticity_modulus: 1.0e+09 - # poissons_ratio: 0.3 - # body_force_x: 0 - # body_force_y: 0 - # bc_list: ["left"] diff --git a/input/panels3d.yaml b/input/panels3d.yaml deleted file mode 100644 index a4b5595..0000000 --- a/input/panels3d.yaml +++ /dev/null @@ -1,59 +0,0 @@ -general: - test: False - geometry_module: panels3d - output_dir: output/panels3d - mesh_only: False - # input_mesh_dir: output/panels3d/mesh - structural_analysis: False - fluid_analysis: True -domain: - x_min: -20 - x_max: 100 - y_min: -30 - y_max: 30 # 20+39 39 is panel to panel - z_min: 0 - z_max: 20 - l_char: 4 -pv_array: - stream_rows: 1 - elevation: 1.5 - stream_spacing: 7.0 - panel_chord: 2.0 - panel_span: 7.0 - panel_thickness: 0.1 - tracker_angle: 0 - span_spacing: 12.0 - span_rows: 1 - span_fixation_pts: 3.5 -solver: - dt: .1 #0.0025 - t_final: 10.0 # 10.0 - solver1_ksp: cg - solver2_ksp: cg - solver3_ksp: cg - solver5_ksp: gmres - solver1_pc: jacobi - solver2_pc: jacobi - solver3_pc: jacobi - solver5_pc: jacobi - save_text_interval: .1 #0.01 - save_xdmf_interval: .1 #0.01 -fluid: - u_ref: 1.0 - nu: 0.01 - turbulence_model: smagorinsky #null # - periodic: false - bc_y_max: slip # slip noslip free - bc_y_min: slip # slip noslip free - bc_z_max: slip # slip noslip free - bc_z_min: noslip # slip noslip free - wind_direction: 250 # slip noslip free -structure: - dt : 0.1 - elasticity_modulus: 1.0e+05 - poissons_ratio: 0.3 - body_force_x: 0 - body_force_y: 0 - body_force_z: -1 #100 - bc_list: [left ] - tube_connection: False \ No newline at end of file diff --git a/input/single_heliostat.yaml b/input/single_heliostat.yaml deleted file mode 100644 index c0bd5bf..0000000 --- a/input/single_heliostat.yaml +++ /dev/null @@ -1,62 +0,0 @@ -general: - geometry_module: heliostats3d - output_dir: output/heliostat - mesh_only: false - structural_analysis: true - fluid_analysis: true -domain: - x_min: -20 - x_max: 40 - y_min: -20 - y_max: 20 # 20+39 39 is panel to panel - z_min: 0 - z_max: 40 - l_char: 4 #1.25 # 1.0 -pv_array: - stream_rows: 1 - elevation: 5.5 - stream_spacing: 10.0 # number from chris ivanov, distance from pier to edge of fixed tilt - panel_chord: 10.0 - panel_span: 11.0 - panel_thickness: 0.1 - tracker_angle: 0 - span_spacing: 1.0 - span_rows: 1 - # span_fixation_pts: [0.6, 6.9, 13.2, 18.3, 22.4] - span_fixation_pts: [5.5] -solver: - dt: .01 #0.0025 - t_final: 20.0 #20.0 # 10.0 - solver1_ksp: gmres - solver2_ksp: gmres - solver3_ksp: gmres - solver5_ksp: gmres - solver1_pc: hypre - solver2_pc: hypre - solver3_pc: hypre - solver5_pc: hypre - save_text_interval: .02 #0.01 - save_xdmf_interval: .02 #0.01 -fluid: - # time_varying_inflow_bc: false # true - u_ref: 2.0 - rho: 1.0 - nu: 1.8e-05 - turbulence_model: smagorinsky #null # - bc_y_max: slip # slip noslip free - bc_y_min: slip # slip noslip free - bc_z_max: slip # slip noslip free - bc_z_min: noslip # slip noslip free - wind_direction: 270 # slip noslip free -structure: - dt : 0.01 # 0.1 - rho: 124.0 - poissons_ratio: 0.3 - elasticity_modulus: 1.0e+09 # 1.0e+03 - body_force_x: 0 - body_force_y: 0 - body_force_z: 0 #100 - bc_list: [] - motor_connection: true - tube_connection: true - beta_relaxation: 0.005 From 498555a3a230dc17060d4bfc4b92379a7e4d4a8b Mon Sep 17 00:00:00 2001 From: arswalid Date: Thu, 26 Jun 2025 16:05:10 -0600 Subject: [PATCH 21/28] black format --- .gitignore | 2 +- pvade/fluid/FlowManager.py | 17 +++++++---------- 2 files changed, 8 insertions(+), 11 deletions(-) diff --git a/.gitignore b/.gitignore index 5b92f10..afb29c3 100644 --- a/.gitignore +++ b/.gitignore @@ -8,4 +8,4 @@ *.h5 docs/_build/ -results/ \ No newline at end of file +results/input/ diff --git a/pvade/fluid/FlowManager.py b/pvade/fluid/FlowManager.py index c771071..cdda51d 100644 --- a/pvade/fluid/FlowManager.py +++ b/pvade/fluid/FlowManager.py @@ -80,7 +80,7 @@ def __init__(self, domain, params): # find hmin in mesh num_cells = domain.fluid.msh.topology.index_map(self.ndim).size_local - + h = dolfinx.cpp.mesh.h(domain.fluid.msh, self.ndim, range(num_cells)) # # This value of hmin is local to the mesh portion owned by the process @@ -94,9 +94,9 @@ def __init__(self, domain, params): print(hmin_local) self.hmin = np.zeros(1) self.hmin = self.comm.allreduce(hmin_local, op=MPI.MIN) - - #------------------------------------------------------------------------------ - + + # ------------------------------------------------------------------------------ + # h = dolfinx.cpp.mesh.h(domain.fluid.msh, self.ndim, range(num_cells)) # # This value of hmin is local to the mesh portion owned by the process @@ -972,7 +972,6 @@ def compute_cfl(self): # self.comm.Allreduce(cfl_max_local, self.cfl_max, op=MPI.MAX) # self.cfl_max = self.cfl_max[0] - # Compute local max only if array has values if self.cfl_vec.vector.array.size > 0: cfl_max_local = np.amax(self.cfl_vec.vector.array) @@ -981,13 +980,11 @@ def compute_cfl(self): # Prepare buffer and reduce across all ranks self.cfl_max = np.zeros(1, dtype=np.float64) - self.comm.Allreduce(np.array(cfl_max_local, dtype=np.float64), self.cfl_max, op=MPI.MAX) + self.comm.Allreduce( + np.array(cfl_max_local, dtype=np.float64), self.cfl_max, op=MPI.MAX + ) self.cfl_max = self.cfl_max[0] - - - - def compute_lift_and_drag(self, params, current_time): self.integrated_force_x = [] From cc11fae55a7e4e7a7cd4c9cbf72547eb9b8d0122 Mon Sep 17 00:00:00 2001 From: arswalid Date: Thu, 26 Jun 2025 16:16:25 -0600 Subject: [PATCH 22/28] conftest --- conftest.py | 14 +++++++------- 1 file changed, 7 insertions(+), 7 deletions(-) diff --git a/conftest.py b/conftest.py index dc94fa2..c473048 100644 --- a/conftest.py +++ b/conftest.py @@ -1,15 +1,11 @@ from pathlib import Path - def pytest_addoption(parser): parser.addoption( - "--input-file", - action="store", - default=None, - help="Run test only for a specific input YAML file", + "--input-file", action="store", default=None, + help="Run test only for a specific input YAML file" ) - def pytest_generate_tests(metafunc): if "input_file" not in metafunc.fixturenames: return @@ -19,8 +15,12 @@ def pytest_generate_tests(metafunc): if input_file_arg: metafunc.parametrize("input_file", [Path(input_file_arg)]) else: - input_dir = Path.cwd() / "input" + # 🔧 Always resolve input/ relative to location of conftest.py (PVade/) + this_dir = Path(__file__).resolve().parent # PVade/ + input_dir = this_dir / "input" + all_files = sorted(input_dir.glob("*.yaml")) if not all_files: raise RuntimeError(f"No input files found in {input_dir}") metafunc.parametrize("input_file", all_files) + From 3928f0d8b30a7eb082b710ee676eee1a319a74e6 Mon Sep 17 00:00:00 2001 From: arswalid Date: Thu, 26 Jun 2025 16:22:47 -0600 Subject: [PATCH 23/28] update location of input --- .github/workflows/test_pvade.yaml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/.github/workflows/test_pvade.yaml b/.github/workflows/test_pvade.yaml index 672a436..b7b2cde 100644 --- a/.github/workflows/test_pvade.yaml +++ b/.github/workflows/test_pvade.yaml @@ -34,7 +34,7 @@ jobs: run: PYTHONPATH=. pytest -sv pvade/tests/ - name: test all inputs shell: bash -l {0} - run: pytest -sv test_all_inputs.py + run: pytest -sv --rootdir=./ test_all_inputs.py # Job 2 of 2 - enforce Black formatting formatting: From fcd4c5673bd385fd6c27109ac6c3c769ebdb6ab3 Mon Sep 17 00:00:00 2001 From: arswalid Date: Thu, 26 Jun 2025 16:28:47 -0600 Subject: [PATCH 24/28] black format --- .github/workflows/test_pvade.yaml | 2 +- .gitignore | 1 + conftest.py | 9 ++-- input/2d_cyld.yaml | 43 +++++++++++++++++++ input/3d_cyld.yaml | 47 ++++++++++++++++++++ input/duramat_case_study.yaml | 61 ++++++++++++++++++++++++++ input/flag2d.yaml | 59 +++++++++++++++++++++++++ input/heated_panels2d.yaml | 71 +++++++++++++++++++++++++++++++ input/inflow_input.yaml | 63 +++++++++++++++++++++++++++ input/panels2d.yaml | 58 +++++++++++++++++++++++++ input/panels3d.yaml | 59 +++++++++++++++++++++++++ input/single_heliostat.yaml | 62 +++++++++++++++++++++++++++ 12 files changed, 531 insertions(+), 4 deletions(-) create mode 100644 input/2d_cyld.yaml create mode 100644 input/3d_cyld.yaml create mode 100644 input/duramat_case_study.yaml create mode 100644 input/flag2d.yaml create mode 100644 input/heated_panels2d.yaml create mode 100644 input/inflow_input.yaml create mode 100644 input/panels2d.yaml create mode 100644 input/panels3d.yaml create mode 100644 input/single_heliostat.yaml diff --git a/.github/workflows/test_pvade.yaml b/.github/workflows/test_pvade.yaml index b7b2cde..897d46c 100644 --- a/.github/workflows/test_pvade.yaml +++ b/.github/workflows/test_pvade.yaml @@ -34,7 +34,7 @@ jobs: run: PYTHONPATH=. pytest -sv pvade/tests/ - name: test all inputs shell: bash -l {0} - run: pytest -sv --rootdir=./ test_all_inputs.py + run: pytest -sv test_all_inputs.py # Job 2 of 2 - enforce Black formatting formatting: diff --git a/.gitignore b/.gitignore index afb29c3..65bed47 100644 --- a/.gitignore +++ b/.gitignore @@ -9,3 +9,4 @@ docs/_build/ results/input/ +input/ diff --git a/conftest.py b/conftest.py index c473048..05a95be 100644 --- a/conftest.py +++ b/conftest.py @@ -1,11 +1,15 @@ from pathlib import Path + def pytest_addoption(parser): parser.addoption( - "--input-file", action="store", default=None, - help="Run test only for a specific input YAML file" + "--input-file", + action="store", + default=None, + help="Run test only for a specific input YAML file", ) + def pytest_generate_tests(metafunc): if "input_file" not in metafunc.fixturenames: return @@ -23,4 +27,3 @@ def pytest_generate_tests(metafunc): if not all_files: raise RuntimeError(f"No input files found in {input_dir}") metafunc.parametrize("input_file", all_files) - diff --git a/input/2d_cyld.yaml b/input/2d_cyld.yaml new file mode 100644 index 0000000..1dc4c57 --- /dev/null +++ b/input/2d_cyld.yaml @@ -0,0 +1,43 @@ +general: + geometry_module: cylinder2d + output_dir: output/cylinder2d + mesh_only: false + structural_analysis: false + fluid_analysis: true +domain: + x_min: 0.0 + x_max: 2.5 + y_min: 0.0 + y_max: 0.41 + l_char: .05 +# pv_array: +# stream_rows: 1 # not used +# elevation: 1.5 # not used +# stream_spacing: 7.0 # not used +# panel_chord: 2.0 # not used +# panel_span: 7.0 # not used +# panel_thickness: 0.1 # not used +# tracker_angle: 30.0 # not used +solver: + dt: 0.000625 + t_final: 3 + solver1_ksp: gmres + solver2_ksp: cg + solver3_ksp: cg + solver1_pc: jacobi + solver2_pc: jacobi + solver3_pc: jacobi + save_text_interval: 0.1 + save_xdmf_interval: 0.1 +fluid: + velocity_profile_type: parabolic + inflow_coeff: 1.5307337294603591 # 4*np.sin(np.pi / 8) # only necessary/used when applying a parabolic vel profile + u_ref: 2.5 + nu: 0.001 # Establish re = 20 with diam = 0.1 and u_bar = u_ref * (4/9) + dpdx: 0.0 + turbulence_model: null + periodic: false + bc_y_max: noslip # slip noslip free + bc_y_min: noslip # slip noslip free + + diff --git a/input/3d_cyld.yaml b/input/3d_cyld.yaml new file mode 100644 index 0000000..8e56c2a --- /dev/null +++ b/input/3d_cyld.yaml @@ -0,0 +1,47 @@ +general: + geometry_module: cylinder3d + output_dir: output/cylinder3d + mesh_only: false + structural_analysis: false + fluid_analysis: true +domain: + x_min: 0.0 + x_max: 2.5 + y_min: 0.0 + y_max: 0.41 + z_min: 0.0 + z_max: 0.41 + l_char: 5.0 #.02 +# pv_array: +# stream_rows: 1 # not used +# elevation: 1.5 # not used +# stream_spacing: 7.0 # not used +# panel_chord: 2.0 # not used +# panel_span: 7.0 # not used +# panel_thickness: 0.1 # not used +# tracker_angle: 30.0 # not used +solver: + dt: 0.001 + t_final: 0.01 + solver1_ksp: gmres + solver2_ksp: cg + solver3_ksp: cg + solver1_pc: jacobi + solver2_pc: jacobi + solver3_pc: jacobi + save_text_interval: 0.001 + save_xdmf_interval: 0.001 +fluid: + velocity_profile_type: parabolic + inflow_coeff: 16 # only necessary/used when applying a parabolic vel profile + u_ref: 0.45 + nu: 0.001 # Establish re = 20 with diam = 0.1 and u_bar = u_ref * (4/9) + dpdx: 0.0 + turbulence_model: null + periodic: false + bc_y_max: noslip # slip noslip free + bc_y_min: noslip # slip noslip free + bc_z_max: noslip # slip noslip free + bc_z_min: noslip # slip noslip free + + diff --git a/input/duramat_case_study.yaml b/input/duramat_case_study.yaml new file mode 100644 index 0000000..459c10b --- /dev/null +++ b/input/duramat_case_study.yaml @@ -0,0 +1,61 @@ +general: + geometry_module: panels3d + output_dir: output/duramat_case_study + mesh_only: false + structural_analysis: true + fluid_analysis: true +domain: + x_min: -20.0 + x_max: 100.0 + y_min: -30.0 + y_max: 30.0 + z_min: 0.0 + z_max: 20.0 + l_char: 1.25 +pv_array: + stream_rows: 1 + stream_spacing: 10.0 + span_rows: 1 + span_spacing: 30.0 + panel_chord: 4.1 + panel_span: 24.25 + panel_thickness: 0.1 + elevation: 2.1 + tracker_angle: 0.0 + span_fixation_pts: [13.2] +solver: + dt: 0.01 + t_final: 20.0 + solver1_ksp: gmres + solver2_ksp: gmres + solver3_ksp: gmres + solver5_ksp: gmres + solver1_pc: hypre + solver2_pc: hypre + solver3_pc: hypre + solver5_pc: hypre + save_text_interval: 0.02 + save_xdmf_interval: 0.02 +fluid: + time_varying_inflow_window: 0.0 + u_ref: 16.0 + rho: 1.0 + nu: 1.8e-05 + turbulence_model: smagorinsky + bc_y_max: slip + bc_y_min: slip + bc_z_max: slip + bc_z_min: noslip + wind_direction: 270.0 +structure: + dt : 0.01 + rho: 124.0 + poissons_ratio: 0.3 + elasticity_modulus: 4.0e+09 + body_force_x: 0.0 + body_force_y: 0.0 + body_force_z: 0.0 + bc_list: [] + motor_connection: true + tube_connection: true + beta_relaxation: 0.5 \ No newline at end of file diff --git a/input/flag2d.yaml b/input/flag2d.yaml new file mode 100644 index 0000000..4312aa8 --- /dev/null +++ b/input/flag2d.yaml @@ -0,0 +1,59 @@ +general: + geometry_module: flag2d + output_dir: output/flag2d + mesh_only: false + structural_analysis: True + fluid_analysis: True +domain: + x_min: 0.0 + x_max: 2.5 + y_min: 0.0 + y_max: 0.41 + l_char: 0.006 # .01 +pv_array: + panel_chord: 0.35 # Sets the length of the flag (starts outside circle radius) + # panel_chord: 0.0 # Sets the length of the flag (starts outside circle radius) + panel_span: 0.05 # Sets the radius of the flag pole + panel_thickness: 0.02 # Sets the thickness of the flag + span_rows: 1 + stream_rows: 1 +solver: + dt: 0.001 # 0.002 + t_final: 30.0 + solver1_ksp: gmres + solver2_ksp: gmres + solver3_ksp: gmres + solver5_ksp: gmres + solver1_pc: hypre + solver2_pc: hypre + solver3_pc: hypre + solver5_pc: hypre + save_text_interval: 0.02 + save_xdmf_interval: 0.02 +fluid: + velocity_profile_type: parabolic + inflow_coeff: 6 # only necessary/used when applying a parabolic vel profile + u_ref: 1.0 # 0.2 1.0 2.0 + time_varying_inflow_window: 2.0 + rho: 1000.0 # 0.2 1.0 2.0 + nu: 0.001 # Establish re = 20 with diam = 0.1 and u = u_ref + dpdx: 0.0 + turbulence_model: null + periodic: false + bc_y_max: noslip # slip noslip free + bc_y_min: noslip # slip noslip free + # warm_up_time: 0.25 # slip noslip free +structure: + dt : 0.001 # 0.002 + # rho : 10000.0 + rho : 10000.0 # 10000.0 + poissons_ratio: 0.4 + elasticity_modulus: 1.4e+06 + body_force_x: 0 + body_force_y: 0 + body_force_z: 0 #100 + bc_list: ["left"] + motor_connection: False + tube_connection: False + beta_relaxation: 0.005 + diff --git a/input/heated_panels2d.yaml b/input/heated_panels2d.yaml new file mode 100644 index 0000000..5c2e85e --- /dev/null +++ b/input/heated_panels2d.yaml @@ -0,0 +1,71 @@ +general: + geometry_module: panels2d + output_dir: output/heatedpanels2d + mesh_only: false + structural_analysis: True + fluid_analysis: True + thermal_analysis: True + debug_flag: True +domain: + x_min: -10 + x_max: 75 + y_min: 0 + y_max: 20 + l_char: 0.1 +pv_array: + stream_rows: 8 + span_rows: 1 + elevation: 1.5 + stream_spacing: 5.0 + panel_chord: 2.0 + panel_span: 7.0 + panel_thickness: 0.1 + tracker_angle: 30. +solver: + dt: .01 + t_final: 60.0 + solver1_ksp: gmres + solver2_ksp: cg + solver3_ksp: cg + solver4_ksp: gmres + solver1_pc: jacobi + solver2_pc: jacobi + solver3_pc: jacobi + solver4_pc: lu + save_text_interval: 0.01 # must be same as or bigger than dt + save_xdmf_interval: 0.01 +fluid: + velocity_profile_type: uniform # loglaw + time_varying_inflow_window: 0.0 + initialize_with_inflow_bc: true + u_ref: 1.0 + nu: 0.001 #15.89e-5 #0.001 + g: 9.81 + beta: 0.00333 + alpha: 2.25e-05 # high Pe # 0.00225 # low Pe (no stability needed) # + turbulence_model: + periodic: false + bc_y_max: slip # slip noslip free + bc_y_min: noslip # slip noslip free + T_ambient: 300.0 + T_bottom: 320.0 + T0_panel: 320.0 +structure: + dt : .01 # 0.002 + # rho : 10000.0 + # rho : 10000.0 # 10000.0 + poissons_ratio: 0.3 + elasticity_modulus: 1.0e+05 + body_force_x: 0 + body_force_y: -1 + body_force_z: 0 #100 + bc_list: ["left"] + motor_connection: False + tube_connection: False + beta_relaxation: 0.005 + + # elasticity_modulus: 1.0e+09 + # poissons_ratio: 0.3 + # body_force_x: 0 + # body_force_y: 0 + # bc_list: ["left"] diff --git a/input/inflow_input.yaml b/input/inflow_input.yaml new file mode 100644 index 0000000..6f3baf8 --- /dev/null +++ b/input/inflow_input.yaml @@ -0,0 +1,63 @@ +general: + geometry_module: flag2d + output_dir: output/inflow + mesh_only: false + structural_analysis: true + fluid_analysis: true +domain: + x_min: 0.0 + x_max: 2.5 + y_min: 0.0 + y_max: 0.41 + l_char: 0.006 # .01 +pv_array: + panel_chord: 0.35 # Sets the length of the flag (starts outside circle radius) + # panel_chord: 0.0 # Sets the length of the flag (starts outside circle radius) + panel_span: 0.01 # Sets the radius of the flag pole + panel_thickness: 0.02 # Sets the thickness of the flag + span_rows: 1 + stream_rows: 1 + elevation: 0.205 +solver: + dt: 0.001 # 0.002 + t_final: 4.0 + solver1_ksp: gmres + solver2_ksp: gmres + solver3_ksp: gmres + # solver4_ksp: gmres + solver1_pc: hypre + solver2_pc: hypre + solver3_pc: hypre + solver4_pc: hypre + save_text_interval: 0.001 # 0.02 + save_xdmf_interval: 0.001 # 0.02 +fluid: + velocity_profile_type: parabolic #loglaw + inflow_coeff: 6 # only necessary/used when applying a parabolic vel profile + u_ref: 4.0 # 0.2 1.0 2.0 + z0: 0.005 #m + d0: 0.05 #m + time_varying_inflow_window: 0.0 + initialize_with_inflow_bc: true + rho: 1000.0 # 0.2 1.0 2.0 + nu: 0.001 # Establish re = 20 with diam = 0.1 and u = u_ref + dpdx: 0.0 + turbulence_model: null + periodic: false + bc_y_max: noslip # slip noslip free + bc_y_min: noslip # slip noslip free + # warm_up_time: 0.25 # slip noslip free +structure: + dt : 0.001 # 0.002 + # rho : 10000.0 + rho : 10000.0 # 10000.0 + poissons_ratio: 0.4 + elasticity_modulus: 1.4e+06 + body_force_x: 0 + body_force_y: 0 + body_force_z: 0 #100 + bc_list: ["left"] + motor_connection: False + tube_connection: False + beta_relaxation: 0.005 + diff --git a/input/panels2d.yaml b/input/panels2d.yaml new file mode 100644 index 0000000..7482b21 --- /dev/null +++ b/input/panels2d.yaml @@ -0,0 +1,58 @@ +general: + geometry_module: panels2d + output_dir: output/panels2d + mesh_only: false + structural_analysis: True + fluid_analysis: True +domain: + x_min: -10 + x_max: 50 + y_min: 0 + y_max: 20 + l_char: 0.5 +pv_array: + stream_rows: 1 + span_rows: 1 + elevation: 1.5 + stream_spacing: 7.0 + panel_chord: 2.0 + panel_span: 7.0 + panel_thickness: 0.1 + tracker_angle: 30. #30.0 +solver: + dt: .001 + t_final: 10.0 + solver1_ksp: gmres + solver2_ksp: cg + solver3_ksp: cg + solver1_pc: jacobi + solver2_pc: jacobi + solver3_pc: jacobi + save_text_interval: 0.001 + save_xdmf_interval: 0.001 +fluid: + u_ref: 0.8 + nu: 0.001 + turbulence_model: + periodic: false + bc_y_max: slip # slip noslip free + bc_y_min: noslip # slip noslip free +structure: + dt : .001 # 0.002 + # rho : 10000.0 + # rho : 10000.0 # 10000.0 + poissons_ratio: 0.3 + elasticity_modulus: 1.0e+05 + body_force_x: 0 + body_force_y: -1 + body_force_z: 0 #100 + bc_list: ["top"] + motor_connection: False + tube_connection: False + beta_relaxation: 0.005 + + # elasticity_modulus: 1.0e+09 + # poissons_ratio: 0.3 + # body_force_x: 0 + # body_force_y: 0 + # bc_list: ["left"] diff --git a/input/panels3d.yaml b/input/panels3d.yaml new file mode 100644 index 0000000..a4b5595 --- /dev/null +++ b/input/panels3d.yaml @@ -0,0 +1,59 @@ +general: + test: False + geometry_module: panels3d + output_dir: output/panels3d + mesh_only: False + # input_mesh_dir: output/panels3d/mesh + structural_analysis: False + fluid_analysis: True +domain: + x_min: -20 + x_max: 100 + y_min: -30 + y_max: 30 # 20+39 39 is panel to panel + z_min: 0 + z_max: 20 + l_char: 4 +pv_array: + stream_rows: 1 + elevation: 1.5 + stream_spacing: 7.0 + panel_chord: 2.0 + panel_span: 7.0 + panel_thickness: 0.1 + tracker_angle: 0 + span_spacing: 12.0 + span_rows: 1 + span_fixation_pts: 3.5 +solver: + dt: .1 #0.0025 + t_final: 10.0 # 10.0 + solver1_ksp: cg + solver2_ksp: cg + solver3_ksp: cg + solver5_ksp: gmres + solver1_pc: jacobi + solver2_pc: jacobi + solver3_pc: jacobi + solver5_pc: jacobi + save_text_interval: .1 #0.01 + save_xdmf_interval: .1 #0.01 +fluid: + u_ref: 1.0 + nu: 0.01 + turbulence_model: smagorinsky #null # + periodic: false + bc_y_max: slip # slip noslip free + bc_y_min: slip # slip noslip free + bc_z_max: slip # slip noslip free + bc_z_min: noslip # slip noslip free + wind_direction: 250 # slip noslip free +structure: + dt : 0.1 + elasticity_modulus: 1.0e+05 + poissons_ratio: 0.3 + body_force_x: 0 + body_force_y: 0 + body_force_z: -1 #100 + bc_list: [left ] + tube_connection: False \ No newline at end of file diff --git a/input/single_heliostat.yaml b/input/single_heliostat.yaml new file mode 100644 index 0000000..c0bd5bf --- /dev/null +++ b/input/single_heliostat.yaml @@ -0,0 +1,62 @@ +general: + geometry_module: heliostats3d + output_dir: output/heliostat + mesh_only: false + structural_analysis: true + fluid_analysis: true +domain: + x_min: -20 + x_max: 40 + y_min: -20 + y_max: 20 # 20+39 39 is panel to panel + z_min: 0 + z_max: 40 + l_char: 4 #1.25 # 1.0 +pv_array: + stream_rows: 1 + elevation: 5.5 + stream_spacing: 10.0 # number from chris ivanov, distance from pier to edge of fixed tilt + panel_chord: 10.0 + panel_span: 11.0 + panel_thickness: 0.1 + tracker_angle: 0 + span_spacing: 1.0 + span_rows: 1 + # span_fixation_pts: [0.6, 6.9, 13.2, 18.3, 22.4] + span_fixation_pts: [5.5] +solver: + dt: .01 #0.0025 + t_final: 20.0 #20.0 # 10.0 + solver1_ksp: gmres + solver2_ksp: gmres + solver3_ksp: gmres + solver5_ksp: gmres + solver1_pc: hypre + solver2_pc: hypre + solver3_pc: hypre + solver5_pc: hypre + save_text_interval: .02 #0.01 + save_xdmf_interval: .02 #0.01 +fluid: + # time_varying_inflow_bc: false # true + u_ref: 2.0 + rho: 1.0 + nu: 1.8e-05 + turbulence_model: smagorinsky #null # + bc_y_max: slip # slip noslip free + bc_y_min: slip # slip noslip free + bc_z_max: slip # slip noslip free + bc_z_min: noslip # slip noslip free + wind_direction: 270 # slip noslip free +structure: + dt : 0.01 # 0.1 + rho: 124.0 + poissons_ratio: 0.3 + elasticity_modulus: 1.0e+09 # 1.0e+03 + body_force_x: 0 + body_force_y: 0 + body_force_z: 0 #100 + bc_list: [] + motor_connection: true + tube_connection: true + beta_relaxation: 0.005 From ec4e351f96db6b8ce6a4e0e653b7132506d3a300 Mon Sep 17 00:00:00 2001 From: arswalid Date: Thu, 26 Jun 2025 16:48:39 -0600 Subject: [PATCH 25/28] reduce from 8 to 4 cores for test all --- test_all_inputs.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test_all_inputs.py b/test_all_inputs.py index 1c0b242..cf603d4 100644 --- a/test_all_inputs.py +++ b/test_all_inputs.py @@ -4,7 +4,7 @@ @pytest.mark.parametrize("mesh_only", [True, False]) -@pytest.mark.parametrize("nprocs", [1, 8]) +@pytest.mark.parametrize("nprocs", [1, 4]) def test_pvade_run(input_file, mesh_only, nprocs): cmd = [ "mpirun", From 0a05747b726d5fa6950236ecd0b27f014620aa93 Mon Sep 17 00:00:00 2001 From: arswalid Date: Mon, 30 Jun 2025 09:28:24 -0600 Subject: [PATCH 26/28] removing 4 cores runs --- test_all_inputs.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test_all_inputs.py b/test_all_inputs.py index cf603d4..3e16def 100644 --- a/test_all_inputs.py +++ b/test_all_inputs.py @@ -4,7 +4,7 @@ @pytest.mark.parametrize("mesh_only", [True, False]) -@pytest.mark.parametrize("nprocs", [1, 4]) +@pytest.mark.parametrize("nprocs", [1]) def test_pvade_run(input_file, mesh_only, nprocs): cmd = [ "mpirun", From f03a4ab1a8b22dc3c67194bd132377dae00f9de9 Mon Sep 17 00:00:00 2001 From: brookeslawski Date: Mon, 30 Jun 2025 12:02:03 -0700 Subject: [PATCH 27/28] fixing errors and removing unnecessary print statements --- pvade/IO/input_schema.yaml | 4 ++-- pvade/fluid/FlowManager.py | 1 - 2 files changed, 2 insertions(+), 3 deletions(-) diff --git a/pvade/IO/input_schema.yaml b/pvade/IO/input_schema.yaml index d3b4d41..1c99362 100644 --- a/pvade/IO/input_schema.yaml +++ b/pvade/IO/input_schema.yaml @@ -119,14 +119,14 @@ properties: properties: stream_rows: default: 0 - minimum: 1 + minimum: 0 maximum: 10 type: "integer" description: "The number of panel rows in the streamwise direction." units: "unitless" span_rows: default: 0 - minimum: 1 + minimum: 0 maximum: 10 type: "integer" description: "The number of panel rows in the spanwise direction." diff --git a/pvade/fluid/FlowManager.py b/pvade/fluid/FlowManager.py index cdda51d..f31df14 100644 --- a/pvade/fluid/FlowManager.py +++ b/pvade/fluid/FlowManager.py @@ -115,7 +115,6 @@ def __init__(self, domain, params): ) if self.rank == 0: - print(f"hmin on fluid = {self.hmin}") print(f"Total num dofs on fluid= {self.num_Q_dofs + self.num_V_dofs}") def build_boundary_conditions(self, domain, params): From c1a91cd49912041585c7fa09300b269c2247d93d Mon Sep 17 00:00:00 2001 From: arswalid Date: Mon, 30 Jun 2025 13:13:14 -0600 Subject: [PATCH 28/28] addition of structure anlysis in panel3d --- input/panels3d.yaml | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/input/panels3d.yaml b/input/panels3d.yaml index a4b5595..8535cc6 100644 --- a/input/panels3d.yaml +++ b/input/panels3d.yaml @@ -4,7 +4,7 @@ general: output_dir: output/panels3d mesh_only: False # input_mesh_dir: output/panels3d/mesh - structural_analysis: False + structural_analysis: True fluid_analysis: True domain: x_min: -20 @@ -56,4 +56,4 @@ structure: body_force_y: 0 body_force_z: -1 #100 bc_list: [left ] - tube_connection: False \ No newline at end of file + tube_connection: False