From 0d8fd6b8e47abfcebfad745019d5ceba5da3aac6 Mon Sep 17 00:00:00 2001
From: Brian Groenke <brian.groenke@awi.de>
Date: Tue, 28 Jan 2025 14:56:51 +0100
Subject: [PATCH] Strip AD types in steady state T init

---
 src/Physics/Heat/heat_init.jl | 11 ++++++++---
 1 file changed, 8 insertions(+), 3 deletions(-)

diff --git a/src/Physics/Heat/heat_init.jl b/src/Physics/Heat/heat_init.jl
index 9e4795ba..2720e1b0 100644
--- a/src/Physics/Heat/heat_init.jl
+++ b/src/Physics/Heat/heat_init.jl
@@ -67,14 +67,19 @@ function steadystate!(sub::SubSurface, state, T0, Qgeo, maxiters::Int, convergen
     z = cells(state.grid)
     ΔT = Inf
     while i < maxiters && ΔT > convergence_thresh
-        T_new = T0 .+ z*Qgeo ./ state.kc
+        T_new = adstrip.(T0 .+ z*Qgeo ./ state.kc)
         ΔT = maximum(abs.(T_new .- state.T))
-        state.T .= T_new
-        # recompute initial condition and diagnostic state
+        state.T .= one(eltype(state.T))*T_new
+        # recompute initial condition and diagnostic state variables
         initialcondition!(sub, state)
         computediagnostic!(sub, state)
         i += 1
     end
+    # compute final temperature profile
+    @. state.T = T0 + z*Qgeo / state.kc
+    # recompute initial condition and diagnostic state variables
+    initialcondition!(sub, state)
+    computediagnostic!(sub, state)
     # return temperature profile
     return state.T
 end
-- 
GitLab