Skip to content

Commit 5bc511a

Browse files
donerancldonerancl
committed
remove constantTVGasDomain
1 parent 1f37bfd commit 5bc511a

File tree

1 file changed

+261
-0
lines changed

1 file changed

+261
-0
lines changed

‎src/Domain.jl

Lines changed: 261 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -102,6 +102,94 @@ function ConstantTPDomain(; phase::E2, initialconds::Dict{X,X2}, constantspecies
102102
end
103103
export ConstantTPDomain
104104

105+
106+
mutable struct ConstantTVGasDomain{N<:AbstractPhase,S<:Integer,W<:Real, W2<:Real, I<:Integer, Q<:AbstractArray} <: AbstractConstantKDomain
107+
phase::N
108+
indexes::Q #assumed to be in ascending order
109+
parameterindexes::Q
110+
constantspeciesinds::Array{S,1}
111+
T::W
112+
V::W
113+
kfs::Array{W,1}
114+
krevs::Array{W,1}
115+
efficiencyinds::Array{I,1}
116+
Gs::Array{W,1}
117+
rxnarray::Array{Int64,2}
118+
mu::W
119+
diffusivity::Array{W,1}
120+
jacobian::Array{W,2}
121+
sensitivity::Bool
122+
alternativepformat::Bool
123+
jacuptodate::MArray{Tuple{1},Bool,1,1}
124+
t::MArray{Tuple{1},W2,1,1}
125+
p::Array{W,1}
126+
thermovariabledict::Dict{String,Int64}
127+
end
128+
129+
function ConstantTVGasDomain(;phase::E2,initialconds::Dict{X,X2},constantspecies::Array{X3,1}=Array{String,1}(),
130+
sparse::Bool=false,sensitivity::Bool=false) where {E<:Real,E2<:AbstractPhase,Q<:AbstractInterface,W<:Real,X,X2,X3}
131+
#set conditions and initialconditions
132+
T = 0.0
133+
V = 0.0
134+
y0 = zeros(length(phase.species)+1)
135+
spnames = [x.name for x in phase.species]
136+
for (key,val) in initialconds
137+
if key == "T"
138+
T = val
139+
elseif key == "V"
140+
V = val
141+
else
142+
ind = findfirst(isequal(key),spnames)
143+
@assert typeof(ind)<: Integer "$key not found in species list: $spnames"
144+
y0[ind] = val
145+
end
146+
end
147+
148+
149+
@assert T != 0.0
150+
@assert V != 0.0
151+
ns = y0[1:end-1]
152+
N = sum(ns)
153+
154+
if length(constantspecies) > 0
155+
spcnames = getfield.(phase.species,:name)
156+
constspcinds = [findfirst(isequal(k),spcnames) for k in constantspecies]
157+
else
158+
constspcinds = Array{Int64,1}()
159+
end
160+
efficiencyinds = [rxn.index for rxn in phase.reactions if typeof(rxn.kinetics)<:AbstractFalloffRate && length(rxn.kinetics.efficiencies) > 0]
161+
Gs = calcgibbs(phase,T)
162+
if :solvent in fieldnames(typeof(phase)) && typeof(phase.solvent) != EmptySolvent
163+
mu = phase.solvent.mu(T)
164+
else
165+
mu = 0.0
166+
end
167+
if phase.diffusionlimited
168+
diffs = [x(T=T,mu=mu,P=P) for x in getfield.(phase.species,:diffusion)]
169+
else
170+
diffs = Array{typeof(T),1}()
171+
end
172+
P = N*R*T/V
173+
C = P/(R*T)
174+
175+
y0[end] = P
176+
kfs,krevs = getkfkrevs(phase,T,V,C,N,ns,Gs,diffs,P,0.0)
177+
kfsp = deepcopy(kfs)
178+
for ind in efficiencyinds
179+
kfsp[ind] = 1.0
180+
end
181+
p = vcat(deepcopy(Gs),kfsp)
182+
if sparse
183+
jacobian=spzeros(typeof(T),length(phase.species),length(phase.species))
184+
else
185+
jacobian=zeros(typeof(T),length(phase.species),length(phase.species))
186+
end
187+
rxnarray = getreactionindices(phase)
188+
return ConstantTPDomain(phase,[1,length(phase.species),length(phase.species)+1],[1,length(phase.species)+length(phase.reactions)],constspcinds,
189+
T,V,kfs,krevs,efficiencyinds,Gs,rxnarray,mu,diffs,jacobian,sensitivity,false,MVector(false),MVector(0.0),p, Dict("P"=>length(phase.species)+1)), y0, p
190+
end
191+
export ConstantTPDomain
192+
105193
struct ConstantVDomain{N<:AbstractPhase,S<:Integer,W<:Real,W2<:Real,I<:Integer,Q<:AbstractArray} <: AbstractVariableKDomain
106194
phase::N
107195
indexes::Q #assumed to be in ascending order
@@ -820,6 +908,7 @@ end
820908

821909
export FragmentBasedConstantTrhoDomain
822910

911+
<<<<<<< HEAD
823912
mutable struct ConstantTLiqFilmDomain{N<:AbstractPhase,S<:Integer,W<:Real,W2<:Real,I<:Integer,Q<:AbstractArray} <: AbstractConstantKDomain
824913
phase::N
825914
indexes::Q #assumed to be in ascending order
@@ -843,6 +932,172 @@ mutable struct ConstantTLiqFilmDomain{N<:AbstractPhase,S<:Integer,W<:Real,W2<:Re
843932
t::MArray{Tuple{1},W2,1,1}
844933
p::Array{W,1}
845934
thermovariabledict::Dict{String,Int64}
935+
=======
936+
937+
@inline function calcthermo(d::ConstantTVGasDomain{W,Y},y::J,t::Q,p::W3=SciMLBase.NullParameters()) where {W3<:SciMLBase.NullParameters,W<:IdealGas,Y<:Integer,J<:Array{Float64,1},Q} #no parameter input
938+
ns = y[d.indexes[1]:d.indexes[2]]
939+
P = y[d.indexes[3]]
940+
N = P*d.V/(R*d.T)
941+
cs = ns./d.V
942+
C = N/d.V
943+
for ind in d.efficiencyinds #efficiency related rates may have changed
944+
d.kfs[ind],d.krevs[ind] = getkfkrev(d.phase.reactions[ind],d.phase,d.T,P,C,N,ns,d.Gs,d.diffusivity,d.V,0.0)
945+
end
946+
return ns,cs,d.T,P,d.V,C,N,d.mu,d.kfs,d.krevs,Array{Float64,1}(),Array{Float64,1}(),Array{Float64,1}(),Array{Float64,1}(),0.0,Array{Float64,1}(),0.0
947+
end
948+
949+
@inline function calcthermo(d::ConstantTVGasDomain{W,Y},y::J,t::Q,p::W2=SciMLBase.NullParameters()) where {W2<:Array{Float64,1},W<:IdealGas,Y<:Integer,J<:Array{Float64,1},Q<:Float64} #uses parameter input
950+
ns = y[d.indexes[1]:d.indexes[2]]
951+
P = y[d.indexes[3]]
952+
N = P*d.V/(R*d.T)
953+
cs = ns./d.V
954+
C = N/d.V
955+
if !d.alternativepformat
956+
@views kfps = p[d.parameterindexes[1]-1+length(d.phase.species)+1:d.parameterindexes[1]-1+length(d.phase.species)+length(d.phase.reactions)]
957+
@views nothermochg= d.Gs == p[d.parameterindexes[1]-1+1:d.parameterindexes[1]-1+length(d.phase.species)]
958+
else
959+
kfps = d.p[length(d.phase.species)+1:end].*p[d.parameterindexes[1]-1+length(d.phase.species)+1:d.parameterindexes[1]-1+length(d.phase.species)+length(d.phase.reactions)]
960+
nothermochg= d.Gs == d.p[1:length(d.phase.species)].+p[d.parameterindexes[1]-1+1:d.parameterindexes[1]-1+length(d.phase.species)]
961+
end
962+
@views nokfchg = count(d.kfs .!= kfps) <= length(d.efficiencyinds) && all(kfps[d.efficiencyinds] .== 1.0)
963+
if nothermochg && nokfchg
964+
for ind in d.efficiencyinds #efficiency related rates may have changed
965+
d.kfs[ind],d.krevs[ind] = getkfkrev(d.phase.reactions[ind],d.phase,d.T,P,C,N,ns,d.Gs,d.diffusivity,d.V,0.0;f=kfps[ind])
966+
end
967+
return ns,cs,d.T,d.P,V,C,N,d.mu,d.kfs,d.krevs,Array{Float64,1}(),Array{Float64,1}(),Array{Float64,1}(),Array{Float64,1}(),0.0,Array{Float64,1}(),0.0
968+
elseif nothermochg
969+
d.kfs = kfps
970+
for ind in d.efficiencyinds #efficiency related rates may have changed
971+
d.kfs[ind],d.krevs[ind] = getkfkrev(d.phase.reactions[ind],d.phase,d.T,P,C,N,ns,d.Gs,d.diffusivity,d.V,0.0;f=kfps[ind])
972+
end
973+
return ns,cs,d.T,P,d.V,C,N,d.mu,d.kfs,d.krevs,Array{Float64,1}(),Array{Float64,1}(),Array{Float64,1}(),Array{Float64,1}(),0.0,Array{Float64,1}(),0.0
974+
else #need to handle thermo changes
975+
d.kfs .= kfps
976+
if !d.alternativepformat
977+
d.Gs = p[d.parameterindexes[1]-1+1:d.parameterindexes[1]-1+length(d.phase.species)]
978+
else
979+
d.Gs = d.p[1:length(d.phase.species)].+p[d.parameterindexes[1]-1+1:d.parameterindexes[1]-1+length(d.phase.species)]
980+
end
981+
krevs = getkfkrevs(d.phase,d.T,P,C,N,ns,d.Gs,d.diffusivity,d.V,0.0;kfs=d.kfs)[2]
982+
for ind in d.efficiencyinds #efficiency related rates may have changed
983+
d.kfs[ind],d.krevs[ind] = getkfkrev(d.phase.reactions[ind],d.phase,d.T,P,C,N,ns,d.Gs,d.diffusivity,d.V,0.0;f=kfps[ind])
984+
end
985+
return ns,cs,d.T,P,d.V,C,N,d.mu,d.kfs,d.krevs,Array{Float64,1}(),Array{Float64,1}(),Array{Float64,1}(),Array{Float64,1}(),0.0,Array{Float64,1}(),0.0
986+
end
987+
end
988+
989+
@inline function calcthermo(d::ConstantTVGasDomain{W,Y},y::Array{W3,1},t::Q,p::W2=SciMLBase.NullParameters()) where {W2,W<:IdealGas,Y<:Integer,W3<:ForwardDiff.Dual,Q} #Autodiff y
990+
ns = y[d.indexes[1]:d.indexes[2]]
991+
P = y[d.indexes[3]]
992+
N = P*d.V/(R*d.T)
993+
cs = ns./d.V
994+
C = N/d.V
995+
if !d.alternativepformat
996+
kfs = convert(typeof(y),p[d.parameterindexes[1]-1+length(d.phase.species)+1:d.parameterindexes[1]-1+length(d.phase.species)+length(d.phase.reactions)])
997+
Gs = p[d.parameterindexes[1]-1+1:d.parameterindexes[1]-1+length(d.phase.species)]
998+
else
999+
kfs = convert(typeof(y),d.p[length(d.phase.species)+1:end].*p[d.parameterindexes[1]-1+length(d.phase.species)+1:d.parameterindexes[1]-1+length(d.phase.species)+length(d.phase.reactions)])
1000+
Gs = d.p[1:length(d.phase.species)].+p[d.parameterindexes[1]-1+1:d.parameterindexes[1]-1+length(d.phase.species)]
1001+
end
1002+
krevs = convert(typeof(y),getkfkrevs(d.phase,d.T,P,C,N,ns,Gs,d.diffusivity,d.V,0.0;kfs=kfs)[2])
1003+
for ind in d.efficiencyinds #efficiency related rates may have changed
1004+
kfs[ind],krevs[ind] = getkfkrev(d.phase.reactions[ind],d.phase,d.T,P,C,N,ns,Gs,d.diffusivity,d.V,0.0;f=kfs[ind])
1005+
end
1006+
return ns,cs,d.T,P,d.V,C,N,d.mu,kfs,krevs,Array{Float64,1}(),Array{Float64,1}(),Array{Float64,1}(),Array{Float64,1}(),0.0,Array{Float64,1}(),0.0
1007+
end
1008+
1009+
@inline function calcthermo(d::ConstantTVGasDomain{W,Y},y::J,t::Q,p::W2=SciMLBase.NullParameters()) where {W2,W<:IdealGas,Y<:Integer,J,Q} #Autodiff p
1010+
ns = y[d.indexes[1]:d.indexes[2]]
1011+
P = y[d.indexes[3]]
1012+
N = P*d.V/(R*d.T)
1013+
cs = ns./d.V
1014+
C = N/d.V
1015+
if !d.alternativepformat
1016+
kfs = p[d.parameterindexes[1]-1+length(d.phase.species)+1:d.parameterindexes[1]-1+length(d.phase.species)+length(d.phase.reactions)]
1017+
Gs = p[d.parameterindexes[1]-1+1:d.parameterindexes[1]-1+length(d.phase.species)]
1018+
else
1019+
kfs = d.p[length(d.phase.species)+1:end].*p[d.parameterindexes[1]-1+length(d.phase.species)+1:d.parameterindexes[1]-1+length(d.phase.species)+length(d.phase.reactions)]
1020+
Gs = d.p[1:length(d.phase.species)].+p[d.parameterindexes[1]-1+1:d.parameterindexes[1]-1+length(d.phase.species)]
1021+
end
1022+
krevs = getkfkrevs(d.phase,d.T,P,C,N,ns,Gs,d.diffusivity,d.V,0.0;kfs=kfs)[2]
1023+
for ind in d.efficiencyinds #efficiency related rates may have changed
1024+
kfs[ind],krevs[ind] = getkfkrev(d.phase.reactions[ind],d.phase,d.T,P,C,N,ns,Gs,d.diffusivity,d.V,0.0;f=kfs[ind])
1025+
end
1026+
return ns,cs,d.T,d.P,V,C,N,d.mu,kfs,krevs,Array{Float64,1}(),Array{Float64,1}(),Array{Float64,1}(),Array{Float64,1}(),0.0,Array{Float64,1}(),0.0
1027+
end
1028+
1029+
@inline function calcthermo(d::ConstantTVGasDomain{W,Y},y::J,t::Q,p::W2=SciMLBase.NullParameters()) where {W2,W<:IdealGas,Y<:Integer,J<:Union{ReverseDiff.TrackedArray,Tracker.TrackedArray},Q} #Autodiff p
1030+
ns = y[d.indexes[1]:d.indexes[2]]
1031+
P = y[d.indexes[3]]
1032+
N = P*d.V/(R*d.T)
1033+
cs = ns./d.V
1034+
C = N/d.V
1035+
kfs = similar(y,length(d.phase.reactions))
1036+
krevs = similar(y,length(d.phase.reactions))
1037+
if !d.alternativepformat
1038+
kfs .= p[d.parameterindexes[1]-1+length(d.phase.species)+1:d.parameterindexes[1]-1+length(d.phase.species)+length(d.phase.reactions)]
1039+
Gs = p[d.parameterindexes[1]-1+1:d.parameterindexes[1]-1+length(d.phase.species)]
1040+
else
1041+
kfs .= d.p[length(d.phase.species)+1:end].*p[d.parameterindexes[1]-1+length(d.phase.species)+1:d.parameterindexes[1]-1+length(d.phase.species)+length(d.phase.reactions)]
1042+
Gs = d.p[1:length(d.phase.species)].+p[d.parameterindexes[1]-1+1:d.parameterindexes[1]-1+length(d.phase.species)]
1043+
end
1044+
krevs .= getkfkrevs(d.phase,d.T,P,C,N,ns,Gs,d.diffusivity,d.V,0.0;kfs=kfs)[2]
1045+
for ind in d.efficiencyinds #efficiency related rates may have changed
1046+
kfs[ind],krevs[ind] = getkfkrev(d.phase.reactions[ind],d.phase,d.T,P,C,N,ns,Gs,d.diffusivity,d.V,0.0;f=kfs[ind])
1047+
end
1048+
return ns,cs,d.T,P,d.V,C,N,d.mu,kfs,krevs,Array{Float64,1}(),Array{Float64,1}(),Array{Float64,1}(),Array{Float64,1}(),0.0,Array{Float64,1}(),0.0
1049+
end
1050+
1051+
@inline function calcthermo(d::ConstantTVGasDomain{W,Y},y::J,t::Q,p::W2=SciMLBase.NullParameters()) where {W2<:Union{ReverseDiff.TrackedArray,Tracker.TrackedArray},W<:IdealGas,Y<:Integer,J,Q} #Tracker/reversediff
1052+
ns = y[d.indexes[1]:d.indexes[2]]
1053+
P = y[d.indexes[3]]
1054+
N = P*d.V/(R*d.T)
1055+
cs = ns./d.V
1056+
C = N/d.V
1057+
if !d.alternativepformat
1058+
Gs = p[d.parameterindexes[1]-1+1:d.parameterindexes[1]-1+length(d.phase.species)]
1059+
kfs = [ind in d.efficiencyinds ? getkfkrev(d.phase.reactions[ind],d.phase,d.T,P,C,N,ns,Gs,d.diffusivity,d.V,0.0)[1]*p[d.parameterindexes[1]-1+length(d.phase.species)+ind] : p[d.parameterindexes[1]-1+length(d.phase.species)+ind] for ind in 1:length(d.phase.reactions)]
1060+
else
1061+
Gs = d.p[1:length(d.phase.species)].+p[d.parameterindexes[1]-1+1:d.parameterindexes[1]-1+length(d.phase.species)]
1062+
kfs = [ind in d.efficiencyinds ? getkfkrev(d.phase.reactions[ind],d.phase,d.T,P,C,N,ns,Gs,d.diffusivity,d.V,0.0)[1]*d.p[length(d.phase.species)+ind]*p[d.parameterindexes[1]-1+length(d.phase.species)+ind] : p[d.parameterindexes[1]-1+length(d.phase.species)+ind] for ind in 1:length(d.phase.reactions)]
1063+
end
1064+
krevs = getkfkrevs(d.phase,d.T,P,C,N,ns,Gs,d.diffusivity,d.V,0.0;kfs=kfs)[2]
1065+
return ns,cs,d.T,P,d.V,C,N,d.mu,kfs,krevs,Array{Float64,1}(),Array{Float64,1}(),Array{Float64,1}(),Array{Float64,1}(),0.0,Array{Float64,1}(),0.0
1066+
end
1067+
1068+
@inline function calcthermo(d::ConstantTVGasDomain{W,Y},y::J,t::Q,p::W2=SciMLBase.NullParameters()) where {W2<:Union{ReverseDiff.TrackedArray,Tracker.TrackedArray},W<:IdealGas,Y<:Integer,J<:Union{ReverseDiff.TrackedArray,Tracker.TrackedArray},Q} #Tracker/reversediff
1069+
ns = y[d.indexes[1]:d.indexes[2]]
1070+
P = y[d.indexes[3]]
1071+
N = P*d.V/(R*d.T)
1072+
cs = ns./d.V
1073+
C = N/d.V
1074+
if !d.alternativepformat
1075+
Gs = p[d.parameterindexes[1]-1+1:d.parameterindexes[1]-1+length(d.phase.species)]
1076+
kfs = [ind in d.efficiencyinds ? getkfkrev(d.phase.reactions[ind],d.phase,d.T,P,C,N,ns,Gs,d.diffusivity,d.V,0.0)[1]*p[d.parameterindexes[1]-1+length(d.phase.species)+ind] : p[d.parameterindexes[1]-1+length(d.phase.species)+ind] for ind in 1:length(d.phase.reactions)]
1077+
else
1078+
Gs = d.p[1:length(d.phase.species)].+p[d.parameterindexes[1]-1+1:d.parameterindexes[1]-1+length(d.phase.species)]
1079+
kfs = [ind in d.efficiencyinds ? getkfkrev(d.phase.reactions[ind],d.phase,d.T,P,C,N,ns,Gs,d.diffusivity,d.V,0.0)[1]*d.p[length(d.phase.species)+ind]*p[d.parameterindexes[1]-1+length(d.phase.species)+ind] : p[d.parameterindexes[1]-1+length(d.phase.species)+ind] for ind in 1:length(d.phase.reactions)]
1080+
end
1081+
krevs = getkfkrevs(d.phase,d.T,P,C,N,ns,Gs,d.diffusivity,d.V,0.0;kfs=kfs)[2]
1082+
return ns,cs,d.T,P,d.V,C,N,d.mu,kfs,krevs,Array{Float64,1}(),Array{Float64,1}(),Array{Float64,1}(),Array{Float64,1}(),0.0,Array{Float64,1}(),0.0
1083+
end
1084+
1085+
1086+
1087+
1088+
1089+
1090+
@inline function calcthermo(d::ConstantTPDomain{W,Y},y::J,t::Q,p::W3=SciMLBase.NullParameters()) where {W3<:SciMLBase.NullParameters,W<:IdealGas,Y<:Integer,J<:Array{Float64,1},Q} #no parameter input
1091+
ns = y[d.indexes[1]:d.indexes[2]]
1092+
V = y[d.indexes[3]]
1093+
N = d.P*V/(R*d.T)
1094+
cs = ns./V
1095+
C = N/V
1096+
for ind in d.efficiencyinds #efficiency related rates may have changed
1097+
d.kfs[ind],d.krevs[ind] = getkfkrev(d.phase.reactions[ind],d.phase,d.T,d.P,C,N,ns,d.Gs,d.diffusivity,V,0.0)
1098+
end
1099+
return ns,cs,d.T,d.P,V,C,N,d.mu,d.kfs,d.krevs,Array{Float64,1}(),Array{Float64,1}(),Array{Float64,1}(),Array{Float64,1}(),0.0,Array{Float64,1}(),0.0
1100+
>>>>>>> 02018d3 (add ConstantTVGasDomain)
8461101
end
8471102

8481103
function ConstantTLiqFilmDomain(; phase::Z, initialconds::Dict{X,E}, constantspecies::Array{X2,1}=Array{String,1}(),
@@ -985,7 +1240,13 @@ end
9851240
return ns, cs, d.T, d.P, V, C, N, d.mu, kfs, krevs, Array{Float64,1}(), Array{Float64,1}(), Array{Float64,1}(), Array{Float64,1}(), 0.0, Array{Float64,1}(), 0.0
9861241
end
9871242

1243+
<<<<<<< HEAD
9881244
@inline function calcthermo(d::ConstantTPDomain{W,Y}, y::J, t::Q, p::W2=SciMLBase.NullParameters()) where {W2,W<:IdealGas,Y<:Integer,J,Q} #Autodiff p
1245+
=======
1246+
1247+
1248+
@inline function calcthermo(d::ConstantTPDomain{W,Y},y::J,t::Q,p::W2=SciMLBase.NullParameters()) where {W2,W<:IdealGas,Y<:Integer,J,Q} #Autodiff p
1249+
>>>>>>> 02018d3 (add ConstantTVGasDomain)
9891250
ns = y[d.indexes[1]:d.indexes[2]]
9901251
V = y[d.indexes[3]]
9911252
N = d.P * V / (R * d.T)

0 commit comments

Comments
 (0)