# Phase field TO in gridap

hi, I am using gridap to do topology optimization along with phase field model. it is mentioned that the adjoint vector should be calculated in each load increament in backward manner. It is confusing for me. by this i should solve the forward solution for displacement and save each output and then again solve for adjoint in backward manner? or i can calculate in each load step while solving for displacement as below

``````delv = 1e-3
const vAppMax = 0.1
innerMax = 10
count = 0
p0 = ones(fem_params.np)
pf_vec = Filter(p0;r, fem_params)
pfh = FEFunction(fem_params.Pf, pf_vec)
pth = (pf -> Threshold(pf; β, η)) ∘ pfh
sPrev = CellState(1.0,dΩ)
sh = project(sPrev ,model,dΩ,order)

w = zeros(2106,730)
wg = zeros(2106,730)

while vApp .< vAppMax
count = count .+ 1
if vApp >= 3e-2
delv = 1e-4
end
vApp = vApp .+ delv
print("\n Entering displacemtent step :", float(vApp))
for  inner = 1: innerMax
ψhPlusPrev = project(ψPlusPrev ,model ,dΩ,order)
RelErr = abs(sum(∫(Gc*ls*∇(sh)⋅∇(sh) + 2*ψhPlusPrev*sh*sh + (Gc/ls)*sh*sh)*dΩ-∫( (Gc/ls)*sh)*dΩ))/abs(sum(∫( (Gc/ls)*sh)*dΩ))
sh_vec = stepPhaseField(fem_params,uh,ψhPlusPrev)
sh = FEFunction(fem_params.U_PF, sh_vec)
u_vec = stepDisp(fem_params,pth, uh, sh, vApp)
uh = FEFunction(fem_params.U0_Disp, u_vec)
ψhPos_in = ψPos∘(ε(uh))
update_state!( new_EnergyState ,ψPlusPrev ,ψhPos_in)
if   RelErr  < 1e-20
break
end

A_mat = MatrixA(pth;fem_params)
G_mat = MatrixG(fem_params)
O_mat = MatrixOf(fem_params)
w_vec =  A_mat' \ (O_mat * u_vec)
w_vecg =  G_mat' \ (Gc * ((sh_vec .- 1)/ls) + ls*sh_vec)
wconjh = FEFunction(fem_params.U0_Disp, conj(w_vec))
wconjhg = FEFunction(fem_params.U_PF, conj(w_vecg))
l_temp(dp) = ∫( (-2*(DAdpf(wconjh,uh , pfh, uh, sh; β, η)))*dp )fem_params.dΩ
l_tempg(dpg) = ∫( -DGdpf(wconjhg,sh , pfh; β, η) *dpg)fem_params.dΩ
dgfdpf = assemble_vector(l_temp, fem_params.Pf)
dggfdpf = assemble_vector(l_tempg, fem_params.Pf)
w[:,count] = dgfdpf
wg[:,count] = dggfdpf
end
end
``````

thank you. i appreciate your suggestions.