diff --git a/src/sipnet/sipnet.c b/src/sipnet/sipnet.c index 17cb6fdbb..9a7d0d8f1 100644 --- a/src/sipnet/sipnet.c +++ b/src/sipnet/sipnet.c @@ -390,13 +390,14 @@ void readParamData(ModelParams **modelParamsPtr, const char *paramFile) { initializeOneModelParam(modelParams, "mineralNInit", &(params.minNInit), ctx.nitrogenCycle); initializeOneModelParam(modelParams, "nVolatilizationFrac", &(params.nVolatilizationFrac), ctx.nitrogenCycle); initializeOneModelParam(modelParams, "nLeachingFrac", &(params.nLeachingFrac), ctx.nitrogenCycle); + initializeOneModelParam(modelParams, "nFixationFrac", &(params.nFixationFrac), ctx.nitrogenCycle); // NOLINTEND - // clang-format on + // clang-format on - readModelParams(modelParams, paramF); + readModelParams(modelParams, paramF); - fclose(paramF); + fclose(paramF); } /*! @@ -410,9 +411,9 @@ void outputHeader(FILE *out) { fprintf(out, "year day time plantWoodC plantLeafC woodCreation "); fprintf(out, "soil microbeC coarseRootC fineRootC "); fprintf(out, "litter soilWater soilWetnessFrac snow "); - fprintf(out, - "npp nee cumNEE gpp rAboveground rSoil rRoot ra rh rtot " - "evapotranspiration fluxestranspiration minN n2oFlux nLeachFlux\n"); + fprintf(out, "npp nee cumNEE gpp rAboveground rSoil rRoot ra rh rtot " + "evapotranspiration fluxestranspiration minN n2oFlux nLeachFlux " + "nFixationFlux\n"); } /*! @@ -433,12 +434,12 @@ void outputState(FILE *out, int year, int day, double time) { trackers.soilWetnessFrac, envi.snow); fprintf(out, "%8.2f %8.2f %8.2f %8.2f %8.3f %8.3f %8.3f %8.3f %8.3f %8.3f %8.8f " - "%8.4f %8.3f %8.6f %8.4f\n", + "%8.4f %8.3f %8.6f %8.4f %8.6f\n", trackers.npp, trackers.nee, trackers.totNee, trackers.gpp, trackers.rAboveground, trackers.rSoil, trackers.rRoot, trackers.ra, trackers.rh, trackers.rtot, trackers.evapotranspiration, fluxes.transpiration, envi.minN, fluxes.nVolatilization, - fluxes.nLeaching); + fluxes.nLeaching, fluxes.nFixation); } // de-allocate space used for climate linked list @@ -1241,6 +1242,16 @@ void calcNLeachingFlux() { fluxes.nLeaching = envi.minN * phi * params.nLeachingFrac; } +/*! + * Calculate mineral N fixation flux + */ +void calcNFixationFlux() { + // flux = k_fix * NPP, g N * m^-2 * day^-1 + // k_fix represents C cost to fix N in range 0.07-0.17 gN/gC + + fluxes.nFixation = params.nFixationFrac * getMeanTrackerMean(meanNPP); +} + /*! * Calculate flux terms for sipnet as part of main model flow * @@ -1310,6 +1321,7 @@ void calculateFluxes(void) { // Leaching depends on drainage flux so makes sure calcNLeachingFlux // occurs after calcSoilWaterFluxes calcNLeachingFlux(); + calcNFixationFlux(); } } @@ -1598,7 +1610,8 @@ void updatePoolsForSoil(void) { climate->length; // Nitrogen Cycle - envi.minN -= (fluxes.nVolatilization + fluxes.nLeaching) * climate->length; + envi.minN += (fluxes.nFixation - fluxes.nVolatilization - fluxes.nLeaching) * + climate->length; } // !!! main runner function !!! diff --git a/src/sipnet/state.h b/src/sipnet/state.h index d107033d7..328e11f1e 100644 --- a/src/sipnet/state.h +++ b/src/sipnet/state.h @@ -367,6 +367,9 @@ typedef struct Parameters { // Fraction of mineral N lost to leaching per day double nLeachingFrac; + // Fraction of NPP available to be fixed by vegetation + double nFixationFrac; + } Params; #define NUM_PARAMS (sizeof(Params) / sizeof(double)) @@ -535,6 +538,8 @@ typedef struct FluxVars { double nVolatilization; // Mineral N lost to leaching double nLeaching; + // Mineral N gained through plant fixation + double nFixation; // **************************************** // Fluxes for event handling diff --git a/tests/sipnet/test_modeling/testNitrogenCycle.c b/tests/sipnet/test_modeling/testNitrogenCycle.c index e3999dc0a..2f4c12df4 100644 --- a/tests/sipnet/test_modeling/testNitrogenCycle.c +++ b/tests/sipnet/test_modeling/testNitrogenCycle.c @@ -17,7 +17,7 @@ void setupTests(void) { ctx.litterPool = 1; ctx.nitrogenCycle = 1; } - +/* void initState(double initN, double nVol, double nLeachFrac) { // per test envi.minN = initN; @@ -40,7 +40,83 @@ void initState(double initN, double nVol, double nLeachFrac) { params.baseSoilResp = 0.06; params.soilRespQ10 = 2.9; } +*/ +void initNLeachingState(double initN, double nLeachFrac) { + // per test + envi.minN = initN; + params.nLeachingFrac = nLeachFrac; + + // static + envi.soilWater = 5.0; + envi.soil = 1.5; + envi.litter = 1; + + params.minNInit = 0.0; + params.soilWHC = 10.0; + + fluxes.nVolatilization = 0; + fluxes.nLeaching = 0; + + // Values from russell_2 smoke test + params.soilRespMoistEffect = 1.0; + params.baseSoilResp = 0.06; + params.soilRespQ10 = 2.9; +} + +void initNVolatilizationState(double initN, double nVol) { + // per test + envi.minN = initN; + params.nVolatilizationFrac = nVol; + + // static + envi.soilWater = 5.0; + envi.soil = 1.5; + envi.litter = 1; + + params.minNInit = 0.0; + params.soilWHC = 10.0; + + fluxes.nVolatilization = 0; + + // Values from russell_2 smoke test + params.soilRespMoistEffect = 1.0; + params.baseSoilResp = 0.06; + params.soilRespQ10 = 2.9; +} + +void initNFixationState(double initN, double nFixFrac) { + // per test + envi.minN = initN; + params.nFixationFrac = nFixFrac; + + // static + envi.soilWater = 5.0; + envi.soil = 1.5; + envi.litter = 1; + + params.minNInit = 0.0; + params.soilWHC = 10.0; + + fluxes.nVolatilization = 0; + fluxes.nLeaching = 0; + fluxes.nFixation = 0; + + // Values from russell_2 smoke test + params.soilRespMoistEffect = 1.0; + params.baseSoilResp = 0.06; + params.soilRespQ10 = 2.9; +} + +int checkFlux(double calcFlux, double expFlux) { + int status = 0; + if (!compareDoubles(calcFlux, expFlux)) { + status = 1; + logTest("Calculated N flux is %f, expected %f\n", calcFlux, expFlux); + } + return status; +} +/* int checkVolatilizationFlux(double expNVolFlux) { int status = 0; if (!compareDoubles(fluxes.nVolatilization, expNVolFlux)) { @@ -62,6 +138,41 @@ int checkNLeachingFlux(double expNLeachingFlux) { return status; } +*/ + +int testNFixation(void) { + int status = 0; + double minN; + double nFixFrac; + double expNFixation; + + minN = 1; + nFixFrac = 0.25; + initNFixationState(minN, nFixFrac); + expNFixation = nFixFrac; + calcNFixationFlux(); + status |= checkFlux(fluxes.nFixation, expNFixation); + + minN = 1; + nFixFrac = 0.5; + initNFixationState(minN, nFixFrac); + expNFixation = nFixFrac; + calcNFixationFlux(); + status |= checkFlux(fluxes.nFixation, expNFixation); + + // Check minN for the last + updatePoolsForSoil(); + double expMinN = minN + (expNFixation * climate->length); + int minStatus = 0; + if (!compareDoubles(envi.minN, expMinN)) { + minStatus = 1; + logTest("Fixation minN pool is %8.3f, expected %8.3f\n", envi.minN, + expMinN); + } + status |= minStatus; + + return status; +} int testNLeaching(void) { int status = 0; @@ -73,7 +184,7 @@ int testNLeaching(void) { minN = 1; nLeachFrac = 0.5; fluxes.drainage = 5; - initState(minN, 0, nLeachFrac); + initNLeachingState(minN, nLeachFrac); if ((fluxes.drainage / params.soilWHC) < 1) { phi = fluxes.drainage / params.soilWHC; } else { @@ -81,12 +192,13 @@ int testNLeaching(void) { } expNLeaching = minN * phi * nLeachFrac; calcNLeachingFlux(); - status |= checkNLeachingFlux(expNLeaching); + // status |= checkNLeachingFlux(expNLeaching); + status |= checkFlux(fluxes.nLeaching, expNLeaching); minN = 1; nLeachFrac = 0.5; fluxes.drainage = 20; - initState(minN, 0, nLeachFrac); + initNLeachingState(minN, nLeachFrac); if ((fluxes.drainage / params.soilWHC) < 1) { phi = fluxes.drainage / params.soilWHC; } else { @@ -94,7 +206,8 @@ int testNLeaching(void) { } expNLeaching = minN * phi * nLeachFrac; calcNLeachingFlux(); - status |= checkNLeachingFlux(expNLeaching); + // status |= checkNLeachingFlux(expNLeaching); + status |= checkFlux(fluxes.nLeaching, expNLeaching); // Check minN for the last updatePoolsForSoil(); @@ -102,7 +215,8 @@ int testNLeaching(void) { int minStatus = 0; if (!compareDoubles(envi.minN, expMinN)) { minStatus = 1; - logTest("minN pool is %8.3f, expected %8.3f\n", envi.minN, expMinN); + logTest("Leaching minN pool is %8.3f, expected %8.3f\n", envi.minN, + expMinN); } status |= minStatus; @@ -117,19 +231,21 @@ int testNVolatilization(void) { minN = 2; nVolFrac = 0.1; - initState(minN, nVolFrac, 0); + initNVolatilizationState(minN, nVolFrac); double tEffect = calcTempEffect(climate->tsoil); double mEffect = calcMoistEffect(envi.soilWater, params.soilWHC); expNVolFlux = nVolFrac * minN * tEffect * mEffect; calcNVolatilizationFlux(); - status |= checkVolatilizationFlux(expNVolFlux); + // status |= checkVolatilizationFlux(expNVolFlux); + status |= checkFlux(fluxes.nVolatilization, expNVolFlux); // easy proportionality test - doubling params should double output minN *= 2; expNVolFlux *= 2; - initState(minN, nVolFrac, 0); + initNVolatilizationState(minN, nVolFrac); calcNVolatilizationFlux(); - status |= checkVolatilizationFlux(expNVolFlux); + // status |= checkVolatilizationFlux(expNVolFlux); + status |= checkFlux(fluxes.nVolatilization, expNVolFlux); // Check minN for the last updatePoolsForSoil(); @@ -137,7 +253,8 @@ int testNVolatilization(void) { int minStatus = 0; if (!compareDoubles(envi.minN, expMinN)) { minStatus = 1; - logTest("minN pool is %8.3f, expected %8.3f\n", envi.minN, expMinN); + logTest("Volatilization minN pool is %8.3f, expected %8.3f\n", envi.minN, + expMinN); } status |= minStatus; @@ -151,7 +268,7 @@ int testFertilization(void) { double expMinN, expEventMinNFlux, expNVolFlux; // init minN 2, nVol 0.1 - initState(initN, nVolFrac, 0); + initNVolatilizationState(initN, nVolFrac); // fert event: 15 5 10 double fertMinN = 10; @@ -201,6 +318,8 @@ int run(void) { status |= testNLeaching(); + status |= testNFixation(); + return status; } diff --git a/tests/smoke/russell_2/sipnet.param b/tests/smoke/russell_2/sipnet.param index 0c84ae1fc..75387038b 100644 --- a/tests/smoke/russell_2/sipnet.param +++ b/tests/smoke/russell_2/sipnet.param @@ -80,3 +80,4 @@ microbePulseEff 0.45 mineralNInit 1.0 nVolatilizationFrac 0.1 nLeachingFrac 0.25 +nFixationFrac 0.11 \ No newline at end of file