From efc9d6817457b67cbf9313495676ae844f8e2c31 Mon Sep 17 00:00:00 2001 From: "Z J Wegert (Workstation)" Date: Tue, 17 Dec 2024 15:37:45 +1000 Subject: [PATCH] testing zygote --- .../zygote/thermal_compliance_ALM copy.jl | 76 ++++++++++++++++ .../zygote/thermal_compliance_ALM copy_2.jl | 86 +++++++++++++++++++ 2 files changed, 162 insertions(+) create mode 100644 scripts/Embedded/MWEs/zygote/thermal_compliance_ALM copy.jl create mode 100644 scripts/Embedded/MWEs/zygote/thermal_compliance_ALM copy_2.jl diff --git a/scripts/Embedded/MWEs/zygote/thermal_compliance_ALM copy.jl b/scripts/Embedded/MWEs/zygote/thermal_compliance_ALM copy.jl new file mode 100644 index 0000000..f3f2d1a --- /dev/null +++ b/scripts/Embedded/MWEs/zygote/thermal_compliance_ALM copy.jl @@ -0,0 +1,76 @@ +using Gridap, GridapTopOpt + +path="./results/TEST_thermal_compliance_ALM/" +## Parameters +order = 1 +xmax=ymax=1.0 +prop_Γ_N = 0.2 +prop_Γ_D = 0.2 +dom = (0,xmax,0,ymax) +el_size = (200,200) +γ = 0.1 +γ_reinit = 0.5 +max_steps = floor(Int,order*minimum(el_size)/10) +tol = 1/(5*order^2)/minimum(el_size) +κ = 1 +vf = 0.4 +η_coeff = 2 +α_coeff = 4max_steps*γ +iter_mod = 10 +mkpath(path) + +## FE Setup +model = CartesianDiscreteModel(dom,el_size); +el_Δ = get_el_Δ(model) +f_Γ_D(x) = (x[1] ≈ 0.0 && (x[2] <= ymax*prop_Γ_D + eps() || + x[2] >= ymax-ymax*prop_Γ_D - eps())) +f_Γ_N(x) = (x[1] ≈ xmax && ymax/2-ymax*prop_Γ_N/2 - eps() <= x[2] <= + ymax/2+ymax*prop_Γ_N/2 + eps()) +update_labels!(1,model,f_Γ_D,"Gamma_D") +update_labels!(2,model,f_Γ_N,"Gamma_N") + +## Triangulations and measures +Ω = Triangulation(model) +Γ_N = BoundaryTriangulation(model,tags="Gamma_N") +dΩ = Measure(Ω,2*order) +dΓ_N = Measure(Γ_N,2*order) +vol_D = sum(∫(1)dΩ) + +## Spaces +reffe_scalar = ReferenceFE(lagrangian,Float64,order) +V = TestFESpace(model,reffe_scalar;dirichlet_tags=["Gamma_D"]) +U = TrialFESpace(V,0.0) +V_φ = TestFESpace(model,reffe_scalar) +V_reg = TestFESpace(model,reffe_scalar;dirichlet_tags=["Gamma_N"]) +U_reg = TrialFESpace(V_reg,0) + +## Create FE functions +φh = interpolate(initial_lsf(4,0.2),V_φ) + +## Interpolation and weak form +interp = SmoothErsatzMaterialInterpolation(η = η_coeff*maximum(el_Δ)) +I,H,DH,ρ = interp.I,interp.H,interp.DH,interp.ρ + +a(u,v,φ) = ∫((I ∘ φ)*κ*∇(u)⋅∇(v))dΩ +l(v,φ) = ∫(v)dΓ_N + +## Optimisation functionals +J(u,φ) = ∫((I ∘ φ)*κ*∇(u)⋅∇(u))dΩ +dJ(q,u,φ) = ∫(κ*∇(u)⋅∇(u)*q*(DH ∘ φ)*(norm ∘ ∇(φ)))dΩ; +Vol(u,φ) = ∫(((ρ ∘ φ) - vf)/vol_D)dΩ; +dVol(q,u,φ) = ∫(-1/vol_D*q*(DH ∘ φ)*(norm ∘ ∇(φ)))dΩ + +state_map = AffineFEStateMap(a,l,U,V,V_φ,V_φ,φh) + +# Create StateParamIntegrandWithMeasures +_J = GridapTopOpt.StateParamIntegrandWithMeasure(J,state_map) + +function φ_to_j(φ) + u = state_map(φ) + _J(u,φ) +end + +using Zygote +_data = Zygote.gradient(φ_to_j,get_free_dof_values(φh)) + +sum(_data[1]) \ No newline at end of file diff --git a/scripts/Embedded/MWEs/zygote/thermal_compliance_ALM copy_2.jl b/scripts/Embedded/MWEs/zygote/thermal_compliance_ALM copy_2.jl new file mode 100644 index 0000000..000da8d --- /dev/null +++ b/scripts/Embedded/MWEs/zygote/thermal_compliance_ALM copy_2.jl @@ -0,0 +1,86 @@ +using Gridap, GridapTopOpt, GridapDistributed, PartitionedArrays + +parts = (3,3); +ranks = with_debug() do distribute + distribute(LinearIndices((prod(parts),))) +end + +path="./results/TEST_thermal_compliance_ALM2/" +## Parameters +order = 1 +xmax=ymax=1.0 +prop_Γ_N = 0.2 +prop_Γ_D = 0.2 +dom = (0,xmax,0,ymax) +el_size = (200,200) +γ = 0.1 +γ_reinit = 0.5 +max_steps = floor(Int,order*minimum(el_size)/10) +tol = 1/(5*order^2)/minimum(el_size) +κ = 1 +vf = 0.4 +η_coeff = 2 +α_coeff = 4max_steps*γ +iter_mod = 10 +mkpath(path) + +## FE Setup +model = CartesianDiscreteModel(ranks,parts,dom,el_size); +el_Δ = get_el_Δ(model) +f_Γ_D(x) = (x[1] ≈ 0.0 && (x[2] <= ymax*prop_Γ_D + eps() || + x[2] >= ymax-ymax*prop_Γ_D - eps())) +f_Γ_N(x) = (x[1] ≈ xmax && ymax/2-ymax*prop_Γ_N/2 - eps() <= x[2] <= + ymax/2+ymax*prop_Γ_N/2 + eps()) +update_labels!(1,model,f_Γ_D,"Gamma_D") +update_labels!(2,model,f_Γ_N,"Gamma_N") + +## Triangulations and measures +Ω = Triangulation(model) +Γ_N = BoundaryTriangulation(model,tags="Gamma_N") +dΩ = Measure(Ω,2*order) +dΓ_N = Measure(Γ_N,2*order) +vol_D = sum(∫(1)dΩ) + +## Spaces +reffe_scalar = ReferenceFE(lagrangian,Float64,order) +V = TestFESpace(model,reffe_scalar;dirichlet_tags=["Gamma_D"]) +U = TrialFESpace(V,0.0) +V_φ = TestFESpace(model,reffe_scalar) +V_reg = TestFESpace(model,reffe_scalar;dirichlet_tags=["Gamma_N"]) +U_reg = TrialFESpace(V_reg,0) + +## Create FE functions +φh = interpolate(initial_lsf(4,0.2),V_φ) + +## Interpolation and weak form +interp = SmoothErsatzMaterialInterpolation(η = η_coeff*maximum(el_Δ)) +I,H,DH,ρ = interp.I,interp.H,interp.DH,interp.ρ + +a(u,v,φ) = ∫((I ∘ φ)*κ*∇(u)⋅∇(v))dΩ +l(v,φ) = ∫(v)dΓ_N + +## Optimisation functionals +J(u,φ) = ∫((I ∘ φ)*κ*∇(u)⋅∇(u))dΩ +dJ(q,u,φ) = ∫(κ*∇(u)⋅∇(u)*q*(DH ∘ φ)*(norm ∘ ∇(φ)))dΩ; +Vol(u,φ) = ∫(((ρ ∘ φ) - vf)/vol_D)dΩ; +dVol(q,u,φ) = ∫(-1/vol_D*q*(DH ∘ φ)*(norm ∘ ∇(φ)))dΩ + +state_map = AffineFEStateMap(a,l,U,V,V_φ,V_φ,φh) + +# Create StateParamIntegrandWithMeasures +_J = GridapTopOpt.StateParamIntegrandWithMeasure(J,state_map) + +function φ_to_j(φ) + u = state_map(φ) + _J(u,φ) +end + +using Zygote + +state_map(get_free_dof_values(φh)) + +φ_to_j(get_free_dof_values(φh)) + +_data = Zygote.gradient(φ_to_j,get_free_dof_values(φh)) + +sum(_data[1]) \ No newline at end of file