-
Notifications
You must be signed in to change notification settings - Fork 392
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
fix Super high COP in VRFFluidTCtrl model #10752
base: develop
Are you sure you want to change the base?
Conversation
|
this->TotalCoolingCapacity = TotalCondCoolingCapacity * CoolingPLR; | ||
this->TotalHeatingCapacity = TotalCondHeatingCapacity * HeatingPLR; | ||
this->TotalCoolingCapacity = TotalCondCoolingCapacity * CoolingPLR * CyclingRatio; | ||
this->TotalHeatingCapacity = TotalCondHeatingCapacity * HeatingPLR * CyclingRatio; |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
OK, I can see this one being missed because it's multiplied by PLR and looked correct at first glance. If COP came back in line then it's likely that all the electricity reports are correctly adjusting for cycling. Although it's not a given, because things like defrost, cranckcase heater and evap condenser pump power are small compared to compressor power, and are most likely 0 in the test files. I would check those. What else might have been missed? What happens to piping losses during cycling? What about something like compressor speed? should that also be adjusted by CyclingRatio?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Crankcase heater seems to have considered cycling
VRFRTF = min(1.0, (CyclingRatio / PartLoadFraction));
...
this->CrankCaseHeaterPower = this->CCHeaterPower * (1.0 - VRFRTF);
Defrost seems to have its own "runtime fraction", the FractionalDefrostTime
:
this->DefrostPower = DefrostEIRTempModFac * (this->HeatingCapacity / 1.01667) * FractionalDefrostTime;
} else { // Defrost strategy is resistive
this->DefrostPower = this->DefrostCapacity * FractionalDefrostTime;
}
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Crankcase heater power looks good. I think defrost power needs a CyclingRatio term. FractionalDefrostTime is the defrost time for a full time step. If the system doesn't run that time step (i.e., CyclingRatio = 0.000000001) then DefrostPower should be 0. You can check FractionalDefrostTime and see that it's a function of OA humidity ratio and has nothing to do with PLR or CyclingRatio. Of course you have to check my thinking before doing anything.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I think you're right. FractionalDefrostTime is either from the user input or from OutdoorCoildw
. It should probably multiply by cycling ratio.
For this->EvapCondPumpElecPower
term, since the FluidTCtrl model doesn't have the "Evaporative Condenser Pump Rated Power Consumption" field like the curve-based VRF, this term is always 0.
|
so that when the outdoor conditions cause a non-zero defrost power, AND the system if off, defrost power would be set to zero.
|
|
|
Pdischarge is the compressor discharge pressure, this->CompMaxDeltaP is the maximum compressor pressure rise Pdischarge - this->CompMaxDeltaP is minimum compressor suction pressure RefMinPe is minimum refrigerant evaporating pressure These two terms are both lower bounds of the suction pressure CapMinPe should be no less than both of them, so CapMinPe should be max of the two LB, not min
|
} | ||
// adjust defrost power based on RTF | ||
vrf.DefrostPower *= VRFRTF; |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This looks correct.
@@ -11592,7 +11591,8 @@ void VRFCondenserEquipment::CalcVRFCondenser_FluidTCtrl(EnergyPlusData &state, c | |||
Tdischarge = this->refrig->getSatTemperature(state, max(min(Pdischarge, RefPHigh), RefPLow), RoutineName); | |||
|
|||
// Evaporative capacity ranges_Min | |||
CapMinPe = min(Pdischarge - this->CompMaxDeltaP, RefMinPe); | |||
// suction pressure lower bound need to be no less than both terms in the following | |||
CapMinPe = max(Pdischarge - this->CompMaxDeltaP, RefMinPe); |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This looks correct. Not sure why this didn't show up during the original model develolpment.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Yeah, I'm surprised as well.
@@ -12237,6 +12236,7 @@ void VRFCondenserEquipment::CalcVRFCondenser_FluidTCtrl(EnergyPlusData &state, c | |||
this->ElecHeatingPower = 0; | |||
} | |||
this->VRFCondRTF = VRFRTF; | |||
this->DefrostPower *= VRFRTF; |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This looks correct. A plot of DefrostPower vs RTF or CyclingRatio would show it's working.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This is a plot of heat pump defrost power vs cycling ratio
When the defrost power is non-zero, it is proportional to cycling ratio
The output file is here:
eplusout_defrost_vs_cyclingRatio.xlsx
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Why is defrost electricity so small? at 6 E -9
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I used this test file. Its "Resistive Defrost Heater Capacity" is 0.0000001
US+SF+CZ5B+hp+slab+IECC_2021_VRFPhysics_v2_hrdsize_V2420.expidf.zip
This branch looks good now. I looked at fix_eplusout.err.zip and see this warning. Why is OutputReportTabular looking to convert a cooling coil type to IP units? This is probably a FluidTCtrl cooling coil. You could fix that here or file an issue.
|
@yujiex these warnings (develop_eplusout.err.zip) are set up as recurring warnings so that the err file does not contain thousands of messages. However, this file does contain thousands of messages (300,000 warnings). Can you please investigate, either here or a different branch.
|
Thanks for noticing this. I will fix it here too |
Good point. I will make it recurring warnings |
|
the error message counter is a sum of all previous counter values it overflows to negative non-recurring warning keeps showing up after the overflow use the counter, not the sum of counter in the predicate
@rraustad I've investigated the two issues SI to IP unit conversion warning I fixed the warning of "unit conversion from Cooling Coil Type [1] into IP units", by changing the [] to () in the column header "Cooling Coil Type [1]". The warning is produced because the code assumes there is a unit of measure inside the bracket and will attempt to perform unit conversion. Here this is a categorical variable and should not attempt to do unit conversion, so [] is not appropriate. It is changed to () per suggestion from @JasonGlazer recurring warning shows too many times I also figured out the reason for the recurring warning message not showning as recurring. The messages keeps showing up is not produced by the I believe changing the condition to just
|
|
|
|
|
The fix for the high COP part (#10751) might have some problem. Let me investigate some more. Please hold off on reviewing. Thanks. |
because the demand is roughly Q_evap_req = TU_load + Pipe_Q - Ncomp (multiplied by an adjustment factor C_cap_operation). Ncomp is an input-ouptut variable of the function. Previously Cycling ratio is multiplied outside of VRFOU_CalcCompH after the function has computed Ncomp. However, this will lead to very small cycling ratio non-compatible with the demand as the demand is underestimated because of the Ncomp (in reality it should not be this large, it should be Ncomp * CyclingRatio). This doesn't matter on the cooling side as Ncomp is not involved in the demand calculation there.
|
remove C_cap_operation in demand calculation
|
@yujiex it has been 16 days since this pull request was last updated. |
@yujiex it has been 7 days since this pull request was last updated. |
Pull request overview
1. Fix of the super high COP problem
The following shows the distribution of heating and cooling COP before vs after the fix. For display purpose, the COP=0 points are removed.
eplusout_defect_output_1025.csv
eplusout_fix_1025.csv
COP is calculated as
this->TotalHeatingCapacity
divided bythis->ElecHeatingPower
and other auxiliary power.this->ElecHeatingPower
is calculated as compressor power and outdoor fan powerthis->ElecHeatingPower = this->Ncomp + this->OUFanPower;
The following shows how
this->TotalHeatingCapacity
andthis->TotalCoolingCapacity
are calculated in FluidTCtrl model and the curve-based model. This has led me to believe the super high COP is cause by numerator not multiplying theCyclingRatio
term. But it's NOT because of this.heating COP fix
Even if here the cycling ratio is not multiplied, when plotting the "VRF Heat Pump Total Heating(Cooling) Rate" against cycling ratio, there's actually a positive linear relationship. However, the linear relationship between heating-rate and cycling ratio has a noticeable intercept, while the intercept for the compressor power and cycling ratio line is near zero. This would mean, when cycling ratio is very small, numerator in heating COP will be non-zero while the denominator approaches zero (compressor power and fan power both go to zero), the heating COP will thus blow up.
Comparing with v23.1, this non-zero intercept is also present in the heat pump heating rate and cycling ratio relationship.
However, it didn't blow up as the outdoor unit fan power is always constant and not cycling with the coil even if the "Supply Air Fan Operating Mode Schedule Name" points to a constant-zero schedule (fan always cycles with coil). The following shows the comparison of heatingCOP vs cycling ratio relationship between current develop and v23-1. We can see when zoomed in, they both have the exponential decay pattern.
To fix it, the Ncomp multiplying by cycling ratio is moved to inside of the
VRFOU_CalcCompH
function. Previously (in the develop branch), Ncomp is calculated inside the function by multiplying the rated capacity and the power performance curve value, and the cycling ratio is multiplied to the calculated Ncomp after theVRFOU_CalcCompH
. However, since the heating demand inVRFOU_CalcCompH
involves the compressor heat release (Q_evap_req = TU_load + Pipe_Q - Ncomp;
), if hereNcomp
didn't multiply by cycling ratio, than the compressor heat release is overestimated, leading to an underestimation of the heating demandQ_evap_req
, as a result, the computed cycling ratio will also be smaller than it should be, creating a mismatch between the cycling ratio and the heating rate. Whether the cycling ratio is multiplied inside or afterVRFOU_CalcCompH
will matter in the cycling ratio calculation in heating mode, but not in cooling mode as Ncomp was not in the demand calculation there.After applied this fix, this is how the heating rate, compressor power, and COP look relative to cycling ratio
cooling COP fix
The cooling COP didn't go up as high as heating one but still with a maximum of 150+, which is too high to be realistic. The reason for the high cooling COP is likely the non-linearity between heat pump cooling rate and the cycling ratio when cycling ratio is very small.
In the cycling ratio calculation, the code attempts to match the demand and capacity. If the they cannot match because the capacity is too large, then it will cycle and the cycling ratio is calculated as
When the cycling ratio is very small, this adjustment term
C_cap_operation
exhibits a exponential decay pattern as well and quickly goes to 0 as cycling ratio drops, which seems to be the cause of the non-linear pattern in the heat pump cooling rate-cycling ratio relationship (the following plot is in the develop branch, using debug prints to export value of cyclingRatio and C_cap_operation after line 11486 (Q_h_OU *= CyclingRatio;
line)). The 6 subfigures on the right show how the 6 input variables vary with cycling ratio.h_comp_in_real
,P_evap_reql
, andT_comp_in_rate
show some shoot up very high for when the cycling ratio is low.after removing this multiplier, the cooling rate non-linear behavior is gone and COP vs cycling ratio plot is like the following
Note that in version 23.1, the cycling ratio is like the following
2. Fix of the refrigerant warning
The issue is originated from the change of this line in PR#10416 the following line. The intention is to
The intention of the change is to impose a lower bound in solving for an evaporating temperature when the heating load is small. Such a lower bound is necessary as otherwise, in many time steps when the system is cycling, the OU evaporating temperature will be -72C, which is too low to be feasible. However, using the value of the input field "Variable Evaporating Temperature Minimum for Indoor Unit" as the lower bound will trigger a bunch of refrigerant warnings.
To fix the refrigerant flow problem, it is observed that the following line might have been a mistake
CapMinPe = min(Pdischarge - this->CompMaxDeltaP, RefMinPe);
.Pdischarge
is the compressor discharge pressure,this->CompMaxDeltaP
is the maximum compressor pressure rise. ThenPdischarge - this->CompMaxDeltaP
is minimum compressor suction pressure considering pressure rise bound. TheRefMinPe
is minimum refrigerant evaporating pressure, which is a pressure lower bound considering physical feasibility.CapMinPe
should be no less than both of them, so it should be max of them, not min.In the code base, there are three places where
CapMinPe
is calculated (shown in the following). Two of them are max, one is min, hinting that the min might be a bug.After imposing the right pressure bound, the refrigerant warnings are mostly gone. The following are the eplusout.err files from the develop and the feature branch. A difference can be seen just from the file size.
fix_eplusout.err.zip
develop_eplusout.err.zip
3. Fix of the recurring warning issue
The warning message keep showing up is the non-recurring warning (which should only be shown once or a few times).
The excessive amount of warning is produced by this chunk of code. The condition
this->errors[(int)RefrigError::SatSupDensity].count <= df->RefrigErrorLimitTest
decides whether it prints this warning. The right hand side of the inequality is 1, while the left hand side overflows to negative at some point to -2147450880. From then on, the condition just keeps evaluates to true and the warning keeps showing up. It overflows as it is not a counter, but a sum of all previous counter values. This makes it get large very quickly. When it overflows, the counter value ofdf->SatErrCountGetSupHeatDensityRefrig
is 65536. Thisthis->errors[(int)RefrigError::SatSupDensity].count <= df->RefrigErrorLimitTest
is 1 + 2 + ... + 65536, which overflowed.4. Cooling Coil Type [1] IP to SI unit conversion warning
The warning is produced because the code assumes there is a unit of measure inside the bracket and will attempt to perform unit conversion. Here this is a categorical variable and should not attempt to do unit conversion, so [] is not appropriate. The problem if fixed by changing the "[]" into "()" in this "Cooling Coil Type [1]".
Regression diffs
Table big diff and Table string diffs
There are 733 files have "ColumnHeadingDifference", leading to table big diffs and string diffs. This is caused by changing "Cooling Coil Type [1]" to "Cooling Coil Type (1)" in this feature branch. The diffs are expected.
meter diff
Meter diff happens in test idf "VariableRefrigerantFlow_FluidTCtrl_5Zone", in the following two fields.
The meter diff is caused by this change. The changes here will cause differences in compressor power, which will lead to changes in Electricity:HVAC and Electricity:Facility.
eso diff
in US+SF+CZ4A+hp+crawlspace+IECC_2006_VRF, VariableRefrigerantFlow_FluidTCtrl_5Zone, and VariableRefrigerantFlow_FluidTCtrl_HR_5Zone, the following fields have eso diffs
These diffs are expected as the numerator of the COP term (heating/cooling output) is changed in this feature: it is multiplied by the cycling ratio now. The changes in COP and heating cooling rate are both expected. Changes in Outdoor Unit Evaporating Temperature is due to the correction in refrigerant pressure lower bound, which is used in calculation of OU evaporating temperature. Heating electricity rate and energy difference happen in 1/21. Heating electricity rate is the sum of compressor power and outdoor unit fan power. As a result of the changing
MinOutdoorUnitTe
argument frommax(this->IUEvapTempLow, CapMinTe)
toCapMinTe
in functionVRFOU_CalcCompH
, compressor power will change. So diffs in heating electricity rate is also expected.err diffs
11 files have error diffs. This is because the following warning is removed as a result of the column header change of "Cooling Coil Type".
NOTE: ENHANCEMENTS MUST FOLLOW A SUBMISSION PROCESS INCLUDING A FEATURE PROPOSAL AND DESIGN DOCUMENT PRIOR TO SUBMITTING CODE
Pull Request Author
Add to this list or remove from it as applicable. This is a simple templated set of guidelines.
Reviewer
This will not be exhaustively relevant to every PR.