Seeking advice on model setup

[Rumple]

  Greetings to the esteemed dlparkhurst and all forum users. I am currently attempting to adapt the kinetic model of the Activated Sludge Model ASM1 into a format compatible with PHREEQC(Amm.dat). I would like to seek your advice on the following three questions: 
  How can I ensure that the dissolved oxygen concentration remains at 1.25 mmol/L throughout the KINETICS reaction? I have tried defining it directly in SOLUTION or using EQUILIBRIUM_PHASES, but the results are not satisfactory.

Code: [Select]

SOLUTION_MASTER_SPECIES
    Organic       Organic          0     CH2O         30.0271
    Snd           Snd              0     N             14.0067
    Xnd           Xnd              0     N             14.0067
    Xi            Xi               0     CH2O              30.0271
    Xp            Xp               0     C6H12NO2              130
    Xs            Xs               0     CH2O              30.0271

SOLUTION_SPECIES
Organic = Organic
    log_k     0
    delta_h   0 kJ
Snd = Snd
    log_k     0
    delta_h   0 kJ
Xnd = Xnd
    log_k     0
    delta_h   0 kJ
Xi = Xi
    log_k     0
    delta_h   0 kJ
Xp = Xp
    log_k     0
    delta_h   0 kJ
Xs = Xs
    log_k     0
    delta_h   0 kJ
2 NO3- + 12 H+ + 10 e- = N2 + 6 H2O
   #-log_k   207.08
    -log_k        0

RATES
#p1 Aerobic growth of heterotrophic biomass
        Denitrification_withO
        -start
10 u = (4*32/176)/(24*3600)
20 Ks = 10/32*1e-3
30 Koh = 0.2/32*1e-3
40 f1 = MOL("Organic")/(Ks + MOL("Organic"))
50 f2 = MOL("O2")/(Koh + MOL("O2"))
60 rate =  u * KIN("heterotrophic") * f1 * f2 *176
70 moles = rate * TIME
80 PUT(rate, 1)
90 if (M + moles) < 0 then moles = -M
100 SAVE moles
        -end


#p2 Anoxic growth of heterotrophic biomass
        Denitrification_withoutO
        -start
10 u = (4*32/176)/(24*3600)
20 Ks = 10/32*1e-3
30 Koh = 0.2/32*1e-3
40 Kno =  0.5/14*1e-3
50 ng = 0.8
60 f1 = MOL("Organic")/(Ks + MOL("Organic"))
70 f2 = Koh/(Koh + MOL("O2"))
80 f3 =  TOT("N")/(Kno+TOT("N"))
90 rate = u * KIN("heterotrophic") * f1 * f2 * f3 *ng *176
100 moles = rate * TIME
110 PUT(rate, 2)   # save the rate for use in Biomass rate calculation
120 SAVE moles
        -end


#p3 Aerobic growth of heterotrophic biomass
        Nitrification
        -start
10 u = (0.5*32/176)/(24*3600)
20 Koa = 0.4/32*1e-3
30 Knh = 1/14*1e-3
40 f1 = TOT("Amm")/(Knh + TOT("Amm"))
50 f2 = MOL("O2")/(Koa + MOL("O2"))
60 rate = u * KIN("autotrophic") * f1 * f2 * 176
70 moles = rate * TIME
80 PUT(rate, 3)
90 if (M + moles) < 0 then moles = -M
100 SAVE moles
        -end


#p4  Decay of heterotrophic biomass
        heterotrophic
        -start
10 bh = 0.3/(24*3600)
20 ratewithO = GET(1) 
30 ratewithoutO = GET(2) 
40 rate = -(ratewithO + ratewithoutO)* 0.005681818 + bh*M
50 moles = rate * TIME
60 if (M + moles) < 0 then moles = -M
70 SAVE moles
        -end


#p5  Decay of autotrophic biomass
        autotrophic
        -start
10 b = 0.05/(24*3600)
20 rate = GET(3)
30 rate = -rate*0.005681818 +b*M
40 moles = rate * TIME
50 if (M + moles) < 0 then moles = -M
60 SAVE moles
        -end


#p6 Ammonification
        Ammonification
        -start
10 ka = 0.05*14/(24*3600)
20 rate = ka * MOL("Snd")*KIN("heterotrophic")
30 moles = rate * TIME
40 if (M + moles) < 0 then moles = -M
50 SAVE moles
        -end


#p7 Particulate biodegradable pollution hydrolysis
        organic_slow
        -start
10 kh = 3/(24*3600)
20 Kx = 0.1*(176/32)*1e-3
30 Koh = 0.2/32*1e-3
40 Kno =  0.5/14*1e-3
50 nh = 0.8
60 f1 = (MOL("Xs")*KIN("heterotrophic"))/(Kx+MOL("Xs")/KIN("heterotrophic"))
70 f2 = MOL("O2")/(Koh+ MOL("O2"))
80 f3 = Koh/(Koh+ MOL("O2"))
90 f4 = TOT("N")/(Kno+TOT("N"))
100 rate = kh * f1 * (f2+nh*f3*f4)*KIN("heterotrophic")
110 PUT(rate,4)
120 moles = rate * TIME
130 if (M + moles) < 0 then moles = -M
140 SAVE moles
        -end


#p8 Particulate organic nitrogen hydrolysis
        n_slow
        -start
10 rate = GET(4)
20 rate = rate * MOL("Xnd")/MOL("Xs")
30 moles = rate * TIME
40 if (M + moles) < 0 then moles = -M
50 SAVE moles
        -end
END


SOLUTION 1 Pulse solution with Amm N
    temp      20
    pH        8
    redox     pe
    units     mmol/l
    density   1
    Alkalinity 305 mg/L as HCO3
    Amm       1.285
    Cl        100
    Na        100
    Organic   2.1875
    O(0)      1.25    # O(0)       10    O2(g) 0.3
    Snd       0.357
    Xi        1.25
    Xs        5.3125
    Xnd       0.714
    Xp        1
    -water    1 # kg

#EQUILIBRIUM_PHASE 1
#   O2(g)    -0.3   10

INCREMENTAL_REACTIONS true

KINETICS 1
Denitrification_withO
    -formula  Organic  -0.09375    Amm    -0.005681818   C  0.065340909     H  0.153409091    O  0.082386364
    -tol      1e-08
Denitrification_withoutO
    -formula  Organic  -0.09375    Amm    -0.005681818   N   -0.04995005  C  0.065340909       H    0.153409091     O    0.082386364     Ntg  0.024975025
    -tol      1e-08
Nitrification
    -formula  O   -0.011363636   Amm    -0.303300866      C  -0.028409091    H   0.858766234     N   0.297619048
    -tol      1e-08
autotrophic
    -formula  C5H9NO2  0   Xs    0.02875   Xnd   0.005336991     Xp    0.000344828   C  -0.002409875    H  -0.010501567    O  -0.018076019
    -m        0.000245
    -m0       0.000245
    -tol      1e-08
heterotrophic
    -formula  C5H9NO2  0   Xs    0.02875   Xnd   0.005336991     Xp    0.000344828   C  -0.002409875    H  -0.010501567    O  -0.018076019
    -m        0.004899
    -m0       0.004899
    -tol      1e-08
Ammonification
    -formula  Snd  -0.071428571   Amm  0.071428571
    -tol      1e-08
organic_slow
    -formula        Xs   -0.03125   Organic    0.03125
    -tol      1e-08
n_slow
    -formula        Xnd   -0.071428571   Snd    0.071428571
    -tol      1e-08
-steps       86400 in 24 steps
-step_divide 60
-runge_kutta 3
-bad_step_max 500
END
  In my input, when I only comment out O(0) in SOLUTION 1 (i.e., set the initial concentration of O(g) to 0), I encounter the warning:“ERROR: Numerical method failed on all combinations of convergence parameters, cell/soln/mix 1”. How should I simulate an anaerobic environment?

  In my input, when I only comment out Xp in SOLUTION 1 (i.e., set the initial concentration of Xp to 0), I encounter the warning: "WARNING: Negative moles in solution 1 for Xp, -3.801202e-12. Recovering..." However , from the results when Xp is non-zero, the amount of Xp always increases and should not be consumed. Why does this warning occur? 
 
  I have noticed that many references use units like gCOD/m³ and gN/m³ to represent substance concentrations. Is there a corresponding way to define such units in PHREEQC, or must I manually convert them using conversion factors based on the chemical formula of the substance, as I am currently doing in my input file? 

Code: [Select]

20 Ks = 10/32*1e-3
30 Koh = 0.2/32*1e-3
40 Kno =  0.5/14*1e-3  Any help or suggestions would be greatly appreciated. Thank you!

[dlparkhurst]

The script below maintains the O(0) concentration very close to 1.25 mmol/kgw. Note that concentration is greater than atmospheric oxygen solubility.

WARNINGS are not ERRORS. The warnings appear during the integration calculation, where concentrations are projected to be negative. However, the numerical method backs up and refines the calculation so that negative concentrations are not produced in the final integration. If the numerical methods fail entirely, an ERROR message will be issued and the calculation is aborted.

The calculation works for the oxic environment, you will have to submit an example of what you want to do with an anoxic environment. Note that you have no N in your initial solution, which is would represent NO3-. In an anoxic environment, Amm would be stable, and the pathway from Ntg to Amm tends not to occur.

Code: [Select]

SOLUTION_MASTER_SPECIES
    Organic       Organic          0     CH2O         30.0271
    Snd           Snd              0     N             14.0067
    Xnd           Xnd              0     N             14.0067
    Xi            Xi               0     CH2O              30.0271
    Xp            Xp               0     C6H12NO2              130
    Xs            Xs               0     CH2O              30.0271

SOLUTION_SPECIES
Organic = Organic
    log_k     0
    delta_h   0 kJ
Snd = Snd
    log_k     0
    delta_h   0 kJ
Xnd = Xnd
    log_k     0
    delta_h   0 kJ
Xi = Xi
    log_k     0
    delta_h   0 kJ
Xp = Xp
    log_k     0
    delta_h   0 kJ
Xs = Xs
    log_k     0
    delta_h   0 kJ
2 NO3- + 12 H+ + 10 e- = N2 + 6 H2O
   #-log_k   207.08
    -log_k        0

RATES
#p1 Aerobic growth of heterotrophic biomass
        Denitrification_withO
        -start
10 u = (4*32/176)/(24*3600)
20 Ks = 10/32*1e-3
30 Koh = 0.2/32*1e-3
40 f1 = MOL("Organic")/(Ks + MOL("Organic"))
50 f2 = MOL("O2")/(Koh + MOL("O2"))
60 rate =  u * KIN("heterotrophic") * f1 * f2 *176
70 moles = rate * TIME
80 PUT(rate, 1)
90 if (M + moles) < 0 then moles = -M
100 SAVE moles
        -end


#p2 Anoxic growth of heterotrophic biomass
        Denitrification_withoutO
        -start
10 u = (4*32/176)/(24*3600)
20 Ks = 10/32*1e-3
30 Koh = 0.2/32*1e-3
40 Kno =  0.5/14*1e-3
50 ng = 0.8
60 f1 = MOL("Organic")/(Ks + MOL("Organic"))
70 f2 = Koh/(Koh + MOL("O2"))
80 f3 =  TOT("N")/(Kno+TOT("N"))
90 rate = u * KIN("heterotrophic") * f1 * f2 * f3 *ng *176
100 moles = rate * TIME
110 PUT(rate, 2)   # save the rate for use in Biomass rate calculation
120 SAVE moles
        -end


#p3 Aerobic growth of heterotrophic biomass
        Nitrification
        -start
10 u = (0.5*32/176)/(24*3600)
20 Koa = 0.4/32*1e-3
30 Knh = 1/14*1e-3
40 f1 = TOT("Amm")/(Knh + TOT("Amm"))
50 f2 = MOL("O2")/(Koa + MOL("O2"))
60 rate = u * KIN("autotrophic") * f1 * f2 * 176
70 moles = rate * TIME
80 PUT(rate, 3)
90 if (M + moles) < 0 then moles = -M
100 SAVE moles
        -end


#p4  Decay of heterotrophic biomass
        heterotrophic
        -start
10 bh = 0.3/(24*3600)
20 ratewithO = GET(1)
30 ratewithoutO = GET(2)
40 rate = -(ratewithO + ratewithoutO)* 0.005681818 + bh*M
50 moles = rate * TIME
60 if (M + moles) < 0 then moles = -M
70 SAVE moles
        -end


#p5  Decay of autotrophic biomass
        autotrophic
        -start
10 b = 0.05/(24*3600)
20 rate = GET(3)
30 rate = -rate*0.005681818 +b*M
40 moles = rate * TIME
50 if (M + moles) < 0 then moles = -M
60 SAVE moles
        -end


#p6 Ammonification
        Ammonification
        -start
10 ka = 0.05*14/(24*3600)
20 rate = ka * MOL("Snd")*KIN("heterotrophic")
30 moles = rate * TIME
40 if (M + moles) < 0 then moles = -M
50 SAVE moles
        -end


#p7 Particulate biodegradable pollution hydrolysis
        organic_slow
        -start
10 kh = 3/(24*3600)
20 Kx = 0.1*(176/32)*1e-3
30 Koh = 0.2/32*1e-3
40 Kno =  0.5/14*1e-3
50 nh = 0.8
60 f1 = (MOL("Xs")*KIN("heterotrophic"))/(Kx+MOL("Xs")/KIN("heterotrophic"))
70 f2 = MOL("O2")/(Koh+ MOL("O2"))
80 f3 = Koh/(Koh+ MOL("O2"))
90 f4 = TOT("N")/(Kno+TOT("N"))
100 rate = kh * f1 * (f2+nh*f3*f4)*KIN("heterotrophic")
110 PUT(rate,4)
120 moles = rate * TIME
130 if (M + moles) < 0 then moles = -M
140 SAVE moles
        -end


#p8 Particulate organic nitrogen hydrolysis
        n_slow
        -start
10 rate = GET(4)
20 rate = rate * MOL("Xnd")/MOL("Xs")
30 moles = rate * TIME
40 if (M + moles) < 0 then moles = -M
50 SAVE moles
        -end
END


SOLUTION 1 Pulse solution with Amm N
    temp      20
    pH        8
    redox     pe
    units     mmol/l
    density   1
    Alkalinity 305 mg/L as HCO3
    Amm       1.285
    Cl        100
    Na        100
    Organic   2.1875
    O(0)      1.25   
    Snd       0.357
    Xi        1.25
    Xs        5.3125
    Xnd       0.714
    Xp        1
    -water    1 # kg
END
EQUILIBRIUM_PHASE 1
    O2(g)    -0.34   10
END
INCREMENTAL_REACTIONS true
KINETICS 1
Denitrification_withO
    -formula  Organic  -0.09375    Amm    -0.005681818   C  0.065340909     H  0.153409091    O  0.082386364
    -tol      1e-08
Denitrification_withoutO
    -formula  Organic  -0.09375    Amm    -0.005681818   N   -0.04995005  C  0.065340909       H    0.153409091     O    0.082386364     Ntg  0.024975025
    -tol      1e-08
Nitrification
    -formula  O   -0.011363636   Amm    -0.303300866      C  -0.028409091    H   0.858766234     N   0.297619048
    -tol      1e-08
autotrophic
    -formula  C5H9NO2  0   Xs    0.02875   Xnd   0.005336991     Xp    0.000344828   C  -0.002409875    H  -0.010501567    O  -0.018076019
    -m        0.000245
    -m0       0.000245
    -tol      1e-08
heterotrophic
    -formula  C5H9NO2  0   Xs    0.02875   Xnd   0.005336991     Xp    0.000344828   C  -0.002409875    H  -0.010501567    O  -0.018076019
    -m        0.004899
    -m0       0.004899
    -tol      1e-08
Ammonification
    -formula  Snd  -0.071428571   Amm  0.071428571
    -tol      1e-08
organic_slow
    -formula        Xs   -0.03125   Organic    0.03125
    -tol      1e-08
n_slow
    -formula        Xnd   -0.071428571   Snd    0.071428571
    -tol      1e-08
-steps       86400 in 24 steps
-step_divide 60
-runge_kutta 3
-bad_step_max 500
END
USE solution 1
USE kinetics 1
USE equilibrium_phases 1
END

[Rumple]

 Thank you very much for your guidance. I tried adding N concentration when simulating anoxic conditions, and the calculation ran smoothly. However, If I set Xp to zero, I don't get any error messages, but the calculation keeps running indefinitely without stopping. My initial input should have Xp = 0, but currently, if I set it to zero, the calculation cannot proceed. I need to assign a value to Xp, and different assigned values lead to different results for other substances. What the optimal solution is?
 
 In the -formula definition, for reactions with changes in alkalinity, is it necessary to include expressions like -formula Alkalinity -1, or is it sufficient to just balance the C, H, and O elements in the chemical equation? For instance, in my Ammonification reaction Xnd -> NH4+, which generates 1 mol of alkalinity, which defination is true ?

Code: [Select]

Ammonification
  -formula Snd -1 Amm 1
Ammonification
  -formula Snd -1 AmmH 1
Ammonification
  -formula Snd -1 AmmH 1  OH  1
Ammonification
  -formula Snd -1 AmmH 1  HCO3  1
  Anymore, I don't understand why the PH value has risen when there is more H in the products.
  Simplified input as following.

Code: [Select]

SOLUTION_MASTER_SPECIES
    Organic       Organic          0     CH2O         30.0271
    Snd           Snd              0     N             14.0067
    Xnd           Xnd              0     N             14.0067
    Xi            Xi               0     CH2O              30.0271
    Xp            Xp               0     C6H12NO2              130
    Xs            Xs               0     CH2O              30.0271

SOLUTION_SPECIES
Organic = Organic
    log_k     0
    delta_h   0 kJ
Snd = Snd
    log_k     0
    delta_h   0 kJ
Xnd = Xnd
    log_k     0
    delta_h   0 kJ
Xi = Xi
    log_k     0
    delta_h   0 kJ
Xp = Xp
    log_k     0
    delta_h   0 kJ
Xs = Xs
    log_k     0
    delta_h   0 kJ
2 NO3- + 12 H+ + 10 e- = N2 + 6 H2O
   #-log_k   207.08
    -log_k        0

RATES
#p6 Ammonification
        Ammonification
        -start
10 ka = 0.05*14/(24*3600)
20 rate = ka * MOL("Snd")*KIN("Biomass_yi")
30 moles = rate * TIME
40 if (M + moles) < 0 then moles = -M
50 SAVE moles
        -end
#p4  Decay of heterotrophic biomass
        Biomass_yi
        -start
10 bh = 0.3/(24*3600)
20 ratewithO = GET(1) 
30 ratewithoutO = GET(2) 
40 rate = -(ratewithO + ratewithoutO)* 0.005681818 + bh*M
50 moles = rate * TIME
60 if (M + moles) < 0 then moles = -M
70 SAVE moles
        -end

SOLUTION 1 Pulse solution with Amm N
    temp      20
    pH        8
    redox     pe
    units     mmol/l
    density   1
    Alkalinity 366 mg/L as HCO3
    Amm       0.892857143
    Cl        100
    Na        100
    Organic   2
    O(0)      1.25
    Snd       0.721428571
    Xi        1.25
    Xs        5
    Xnd       1.307142857
    Xp        0.01
    -water    1 # kg

EQUILIBRIUM_PHASES 1
    O2(g)     -0.34   10

SELECTED_OUTPUT
        -file   20260223_2.xlsx
        -high_precision true
USER_PUNCH
        -headings       Organic  Amm    C    H      O  ALK    PH
        -start
  10 PUNCH  TOT("Organic")  TOT("Amm") TOT("C") TOT("H") TOT("O") TOT("Alkalinity")
  20 PUNCH  -LA("H+")
        -end

KINETICS 1
Biomass_yi
    -formula  C5H9NO2  0
    -m        0.004899
    -m0       0.004899
    -tol      1e-08
Ammonification
    -formula  Snd  -0.071428571   Amm   0.071428571
#    -formula  Snd  -0.071428571   AmmH   0.071428571
#    -formula  Snd  -0.071428571   AmmH   0.071428571  OH  0.071428571 
#    -formula  Snd  -0.071428571   AmmH   0.071428571  HCO3  0.071428571 
    -tol      1e-08
-steps       86400 in 24 steps
-step_divide 60
-runge_kutta 3
-bad_step_max 500
END


[dlparkhurst]

You will need to post a better script that demonstrates your problem. The simplified script does not involve Xp. I think in your original, the rates are driving Xp concentration to zero, and I think the consequence is that other reactions cannot proceed. So, my guess is that there is a fundamental flaw in the combination of rates.

I haven't tried to parse through your entire reaction scheme. It looks like Snd is a form of ammonia, so that Snd is converted to a more reactive form of ammonia. In that case, you should define -formula Snd -1 Amm 1. The formula defines the net change in total concentrations of elements. In this case, you are removing Snd from the system and replacing it with Amm. I am assuming there is no change in the concentration of H, C, or O from this reaction, although the speciation of H, C, and O may change due to the reaction.

You should not use Alkalinity anywhere but in a SOLUTION definition.

[Rumple]

Initially, I provided the complete script. Perhaps, as you said, this is an inherent flaw in the combination of these reactions. The only way to address this is for me to modify the rate expressions themselves to prevent this, right?

You are correct. Snd is slowly degradable inorganic nitrogen, and there are changes in alkalinity during this chemical process. However, when I checked the Amm.dat database, I saw that the alkalinity for Amm is 0. The alkalinity I defined for Snd is also 0. If I use the formula `-formula Snd -1 Amm 1`, it does not seem to reflect the change in alkalinity during the chemical reaction. Furthermore, I do not agree with the assumption that the concentrations of H, C, or O remain unchanged in this reaction. If that were the case, the alkalinity should not change either.

[dlparkhurst]

Your original script did not handle oxygen correctly, and may have had other problems.

Sorry, your reaction network is too complicated for me to take the time to sort out. If you can simplify to a simpler part of the system, ammonia release or ammonia oxidation with 2 or 3 kinetic reactions, I will take a look.

I was talking about the single reaction of converting Snd to Amm. That reaction takes place without any change in the C, H, or O concentrations in solution.

In Amm.dat, AmmH+ is defined to have zero alkalinity contribution. Amm, however, will have an alkalinity contribution of 1.0.

Code: [Select]

Amm         AmmH+   0   AmmH            17.031

But lets talk about Snd. Is Snd actually a dissolved species, or does it represent organic matter that degrades to produce ammonia? I am not clear why you need Snd if it is defined similarly to Amm. If n_slow represents the degradation of organic matter, you might define the moles of n_slow (m or KIN("n_slow")) to represent the amount of particulate matter without defining a dissolved species to represent it. The n_slow reaction would simply add Amm (-formula Amm +1) to solution (which would affect alkalinity), and the amount of particulate reactant (M) would decrease. If the moles of kinetic reaction SAVEd is positive, m decreases; if negative, m increases.

Code: [Select]

SOLUTION_MASTER_SPECIES
    Snd           Snd              0     N             14.0067

This defines an alkalinity of zero for Snd, whereas Amm has an alkalinity of 1. So, if you do use -formula Snd -1 Amm +1, then the alkalinity of the solution will increase.

Code: [Select]

SOLUTION_MASTER_SPECIES
    Snd           Snd              0     N             14.0067
SOLUTION_SPECIES
Snd = Snd
    log_k     0
    delta_h   0 kJ

SOLUTION 1 Pulse solution with Amm N
    temp      20
    pH        8
    redox     pe
    units     mmol/l
    density   1
    Alkalinity 305 mg/L as HCO3
    Amm       1.285
    Cl        100
    Na        100
    O(0)      1.25    # O(0)       10    O2(g) 0.3
    Snd       0.357
    -water    1 # kg
REACTION
Snd -1
Amm 1
0.0003 in 5 steps
USER_GRAPH 1
    -axis_titles            "Snd reacted, moles" "Alkalinity" ""
    -initial_solutions      false
    -connect_simulations    true
    -plot_concentration_vs  x
  -start
10 GRAPH_X RXN
20 GRAPH_Y ALK
  -end
    -active                 true
END


[Rumple]

Thank you for your willingness to help. I thought carefully about your previous suggestions and realized that because the autotrophic rate can be positive or negative, it inevitably consumes Xp at the initial stage. When Xp is absent, some warnings appear. Therefore, I attempted to resolve it by splitting the process, and now no warnings occur when Xp = 0.

Code: [Select]

SOLUTION_MASTER_SPECIES
    Organic       Organic          0     CH2O         30.0271
    Snd           Snd              0     N             14.0067
    Xnd           Xnd              0     N             14.0067
    Xi            Xi               0     CH2O              30.0271
    Xp            Xp               0     C6H12NO2              130
    Xs            Xs               0     CH2O              30.0271

SOLUTION_SPECIES
Organic = Organic
    log_k     0
    delta_h   0 kJ
Snd = Snd
    log_k     0
    delta_h   0 kJ
Xnd = Xnd
    log_k     0
    delta_h   0 kJ
Xi = Xi
    log_k     0
    delta_h   0 kJ
Xp = Xp
    log_k     0
    delta_h   0 kJ
Xs = Xs
    log_k     0
    delta_h   0 kJ
2 NO3- + 12 H+ + 10 e- = N2 + 6 H2O
   #-log_k   207.08
    -log_k        0

RATES
#p3 Aerobic growth of heterotrophic biomass
        Nitrification
        -start
10 u = (0.5*32/176)/(24*3600)
20 Koa = 0.4/32*1e-3
30 Knh = 1/14*1e-3
40 f1 = TOT("Amm")/(Knh + TOT("Amm"))
50 f2 = MOL("O2")/(Koa + MOL("O2"))
60 rate = u * KIN("autotrophic") * f1 * f2 * 176
70 moles = rate * TIME
80 PUT(rate, 3)
90 if (M + moles) < 0 then moles = -M
100 SAVE moles
        -end


#p5  Decay of autotrophic biomass
        autotrophic
        -start
10 b = 0.05/(24*3600)
20 rate = GET(3)
30 rate = -rate*0.005681818 +b*M
40 moles = rate * TIME
50 if (M + moles) < 0 then moles = -M
60 SAVE moles
        -end

#p6  Decay of autotrophic biomass
        autotrophic_product
        -start
10 b = 0.05/(24*3600)
20 rate = GET(3)
30 rate = -rate*0.005681818 +b*KIN("autotrophic")
35 if rate < 0 then rate = 0
40 moles = rate * TIME
50 if (M + moles) < 0 then moles = -M
60 SAVE moles
        -end
END

SOLUTION 1 Pulse solution with Amm N
    temp      20
    pH        8
    redox     pe
    units     mmol/l
    density   1
    Alkalinity 305 mg/L as HCO3
    Amm       1.285
    Cl        100
    Na        100
    Organic   2.1875
    O(0)      1.25    # O(0)       10    O2(g) 0.3
    Snd       0.357
    Xi        1.25
    Xs        5.3125
    Xnd       0.714
    # Xp        1
    -water    1 # kg

EQUILIBRIUM_PHASE 1
   O2(g)    -0.3   10

INCREMENTAL_REACTIONS true

KINETICS 1
Nitrification
    -formula  O   -0.011363636   Amm    -0.303300866      C  -0.028409091    H   0.858766234     N   0.297619048
    -tol      1e-08
autotrophic
    -formula  C5H9NO2  0     C  -0.002409875    H  -0.010501567    O  -0.018076019
    -m        0.000245
    -m0       0.000245
    -tol      1e-08
autotrophic_product
    -formula  Xp    0.000344828
    -m        0.000245
    -m0       0.000245
    -tol      1e-08
-steps       86400 in 24 steps
-step_divide 60
-runge_kutta 3
-bad_step_max 500
END
  Sorry, I made a mistake in my previous statement. Snd represents dissolved organic nitrogen, which is produced by the hydrolysis of particulate organic nitrogen and can be converted to ammonia nitrogen through ammonification. Xnd represents particulate biodegradable organic nitrogen, which is generated by microbial decay and removed through the hydrolysis of slowly biodegradable organic matter. I was thinking of the process: Biomass → Xnd → Snd → NH4⁺.

  It turns out that I misunderstood the alkalinity definition of Amm. The testing method you provided was very helpful. Thank you very much.

[dlparkhurst]

You need to conserve charge, so you should not convert Snd (uncharged) to NH4+ (charged) directly.

Snd  (uncharged) ->  Amm  (uncharged) (-formula Snd -1 Amm +1) will avoid some unexpected reactions compared to Snd -> AmmH (-formula Snd -1 AmmH +1), which would essentially be adding Amm + 1/2H2 (a reductant).

[Rumple]

Oh, I see! I think I understand now. Thank you very much.