Fix: account for patch_area in light calculations (#422)

Closes #277 Closes #169 

- ability to set and retrieve patch area
- adjust light calculations based on patch area
- test scm delivers similar results when patch area changes

R/RcppExports.R

@@ -1969,6 +1969,10 @@ Patch___FF16__FF16_Env__size__get <- function(obj_) {
     .Call('_plant_Patch___FF16__FF16_Env__size__get', PACKAGE = 'plant', obj_)
 }
 
+Patch___FF16__FF16_Env__get_area__get <- function(obj_) {
+    .Call('_plant_Patch___FF16__FF16_Env__get_area__get', PACKAGE = 'plant', obj_)
+}
+
 Patch___FF16__FF16_Env__height_max__get <- function(obj_) {
     .Call('_plant_Patch___FF16__FF16_Env__height_max__get', PACKAGE = 'plant', obj_)
 }
@@ -2077,6 +2081,10 @@ Patch___TF24__TF24_Env__size__get <- function(obj_) {
     .Call('_plant_Patch___TF24__TF24_Env__size__get', PACKAGE = 'plant', obj_)
 }
 
+Patch___TF24__TF24_Env__get_area__get <- function(obj_) {
+    .Call('_plant_Patch___TF24__TF24_Env__get_area__get', PACKAGE = 'plant', obj_)
+}
+
 Patch___TF24__TF24_Env__height_max__get <- function(obj_) {
     .Call('_plant_Patch___TF24__TF24_Env__height_max__get', PACKAGE = 'plant', obj_)
 }
@@ -2185,6 +2193,10 @@ Patch___FF16w__FF16_Env__size__get <- function(obj_) {
     .Call('_plant_Patch___FF16w__FF16_Env__size__get', PACKAGE = 'plant', obj_)
 }
 
+Patch___FF16w__FF16_Env__get_area__get <- function(obj_) {
+    .Call('_plant_Patch___FF16w__FF16_Env__get_area__get', PACKAGE = 'plant', obj_)
+}
+
 Patch___FF16w__FF16_Env__height_max__get <- function(obj_) {
     .Call('_plant_Patch___FF16w__FF16_Env__height_max__get', PACKAGE = 'plant', obj_)
 }
@@ -2293,6 +2305,10 @@ Patch___FF16r__FF16_Env__size__get <- function(obj_) {
     .Call('_plant_Patch___FF16r__FF16_Env__size__get', PACKAGE = 'plant', obj_)
 }
 
+Patch___FF16r__FF16_Env__get_area__get <- function(obj_) {
+    .Call('_plant_Patch___FF16r__FF16_Env__get_area__get', PACKAGE = 'plant', obj_)
+}
+
 Patch___FF16r__FF16_Env__height_max__get <- function(obj_) {
     .Call('_plant_Patch___FF16r__FF16_Env__height_max__get', PACKAGE = 'plant', obj_)
 }
@@ -2401,6 +2417,10 @@ Patch___K93__K93_Env__size__get <- function(obj_) {
     .Call('_plant_Patch___K93__K93_Env__size__get', PACKAGE = 'plant', obj_)
 }
 
+Patch___K93__K93_Env__get_area__get <- function(obj_) {
+    .Call('_plant_Patch___K93__K93_Env__get_area__get', PACKAGE = 'plant', obj_)
+}
+
 Patch___K93__K93_Env__height_max__get <- function(obj_) {
     .Call('_plant_Patch___K93__K93_Env__height_max__get', PACKAGE = 'plant', obj_)
 }
@@ -3357,6 +3377,10 @@ StochasticPatch___FF16__FF16_Env__height_max__get <- function(obj_) {
     .Call('_plant_StochasticPatch___FF16__FF16_Env__height_max__get', PACKAGE = 'plant', obj_)
 }
 
+StochasticPatch___FF16__FF16_Env__get_area__get <- function(obj_) {
+    .Call('_plant_StochasticPatch___FF16__FF16_Env__get_area__get', PACKAGE = 'plant', obj_)
+}
+
 StochasticPatch___FF16__FF16_Env__parameters__get <- function(obj_) {
     .Call('_plant_StochasticPatch___FF16__FF16_Env__parameters__get', PACKAGE = 'plant', obj_)
 }
@@ -3441,6 +3465,10 @@ StochasticPatch___TF24__TF24_Env__height_max__get <- function(obj_) {
     .Call('_plant_StochasticPatch___TF24__TF24_Env__height_max__get', PACKAGE = 'plant', obj_)
 }
 
+StochasticPatch___TF24__TF24_Env__get_area__get <- function(obj_) {
+    .Call('_plant_StochasticPatch___TF24__TF24_Env__get_area__get', PACKAGE = 'plant', obj_)
+}
+
 StochasticPatch___TF24__TF24_Env__parameters__get <- function(obj_) {
     .Call('_plant_StochasticPatch___TF24__TF24_Env__parameters__get', PACKAGE = 'plant', obj_)
 }
@@ -3525,6 +3553,10 @@ StochasticPatch___FF16w__FF16_Env__height_max__get <- function(obj_) {
     .Call('_plant_StochasticPatch___FF16w__FF16_Env__height_max__get', PACKAGE = 'plant', obj_)
 }
 
+StochasticPatch___FF16w__FF16_Env__get_area__get <- function(obj_) {
+    .Call('_plant_StochasticPatch___FF16w__FF16_Env__get_area__get', PACKAGE = 'plant', obj_)
+}
+
 StochasticPatch___FF16w__FF16_Env__parameters__get <- function(obj_) {
     .Call('_plant_StochasticPatch___FF16w__FF16_Env__parameters__get', PACKAGE = 'plant', obj_)
 }
@@ -3609,6 +3641,10 @@ StochasticPatch___FF16r__FF16_Env__height_max__get <- function(obj_) {
     .Call('_plant_StochasticPatch___FF16r__FF16_Env__height_max__get', PACKAGE = 'plant', obj_)
 }
 
+StochasticPatch___FF16r__FF16_Env__get_area__get <- function(obj_) {
+    .Call('_plant_StochasticPatch___FF16r__FF16_Env__get_area__get', PACKAGE = 'plant', obj_)
+}
+
 StochasticPatch___FF16r__FF16_Env__parameters__get <- function(obj_) {
     .Call('_plant_StochasticPatch___FF16r__FF16_Env__parameters__get', PACKAGE = 'plant', obj_)
 }
@@ -3693,6 +3729,10 @@ StochasticPatch___K93__K93_Env__height_max__get <- function(obj_) {
     .Call('_plant_StochasticPatch___K93__K93_Env__height_max__get', PACKAGE = 'plant', obj_)
 }
 
+StochasticPatch___K93__K93_Env__get_area__get <- function(obj_) {
+    .Call('_plant_StochasticPatch___K93__K93_Env__get_area__get', PACKAGE = 'plant', obj_)
+}
+
 StochasticPatch___K93__K93_Env__parameters__get <- function(obj_) {
     .Call('_plant_StochasticPatch___K93__K93_Env__parameters__get', PACKAGE = 'plant', obj_)
 }

R/RcppR6.R

@@ -1,6 +1,6 @@
 ## Generated by RcppR6: do not edit by hand
 ## Version: 0.2.4
-## Hash: 1c60048bc1d88826061cc8bb250fbbc8
+## Hash: bbf062e0a16bb3911c8fd115f95cef6b
 
 ##' @importFrom Rcpp evalCpp
 ##' @importFrom R6 R6Class
@@ -2916,6 +2916,13 @@ Patch <- function(T, E) {
           stop("Patch<FF16,FF16_Env>$size is read-only")
         }
       },
+      get_area = function(value) {
+        if (missing(value)) {
+          Patch___FF16__FF16_Env__get_area__get(self)
+        } else {
+          stop("Patch<FF16,FF16_Env>$get_area is read-only")
+        }
+      },
       height_max = function(value) {
         if (missing(value)) {
           Patch___FF16__FF16_Env__height_max__get(self)
@@ -3058,6 +3065,13 @@ Patch <- function(T, E) {
           stop("Patch<TF24,TF24_Env>$size is read-only")
         }
       },
+      get_area = function(value) {
+        if (missing(value)) {
+          Patch___TF24__TF24_Env__get_area__get(self)
+        } else {
+          stop("Patch<TF24,TF24_Env>$get_area is read-only")
+        }
+      },
       height_max = function(value) {
         if (missing(value)) {
           Patch___TF24__TF24_Env__height_max__get(self)
@@ -3200,6 +3214,13 @@ Patch <- function(T, E) {
           stop("Patch<FF16w,FF16_Env>$size is read-only")
         }
       },
+      get_area = function(value) {
+        if (missing(value)) {
+          Patch___FF16w__FF16_Env__get_area__get(self)
+        } else {
+          stop("Patch<FF16w,FF16_Env>$get_area is read-only")
+        }
+      },
       height_max = function(value) {
         if (missing(value)) {
           Patch___FF16w__FF16_Env__height_max__get(self)
@@ -3342,6 +3363,13 @@ Patch <- function(T, E) {
           stop("Patch<FF16r,FF16_Env>$size is read-only")
         }
       },
+      get_area = function(value) {
+        if (missing(value)) {
+          Patch___FF16r__FF16_Env__get_area__get(self)
+        } else {
+          stop("Patch<FF16r,FF16_Env>$get_area is read-only")
+        }
+      },
       height_max = function(value) {
         if (missing(value)) {
           Patch___FF16r__FF16_Env__height_max__get(self)
@@ -3484,6 +3512,13 @@ Patch <- function(T, E) {
           stop("Patch<K93,K93_Env>$size is read-only")
         }
       },
+      get_area = function(value) {
+        if (missing(value)) {
+          Patch___K93__K93_Env__get_area__get(self)
+        } else {
+          stop("Patch<K93,K93_Env>$get_area is read-only")
+        }
+      },
       height_max = function(value) {
         if (missing(value)) {
           Patch___K93__K93_Env__height_max__get(self)
@@ -4823,6 +4858,13 @@ StochasticPatch <- function(T, E) {
           stop("StochasticPatch<FF16,FF16_Env>$height_max is read-only")
         }
       },
+      get_area = function(value) {
+        if (missing(value)) {
+          StochasticPatch___FF16__FF16_Env__get_area__get(self)
+        } else {
+          stop("StochasticPatch<FF16,FF16_Env>$get_area is read-only")
+        }
+      },
       parameters = function(value) {
         if (missing(value)) {
           StochasticPatch___FF16__FF16_Env__parameters__get(self)
@@ -4939,6 +4981,13 @@ StochasticPatch <- function(T, E) {
           stop("StochasticPatch<TF24,TF24_Env>$height_max is read-only")
         }
       },
+      get_area = function(value) {
+        if (missing(value)) {
+          StochasticPatch___TF24__TF24_Env__get_area__get(self)
+        } else {
+          stop("StochasticPatch<TF24,TF24_Env>$get_area is read-only")
+        }
+      },
       parameters = function(value) {
         if (missing(value)) {
           StochasticPatch___TF24__TF24_Env__parameters__get(self)
@@ -5055,6 +5104,13 @@ StochasticPatch <- function(T, E) {
           stop("StochasticPatch<FF16w,FF16_Env>$height_max is read-only")
         }
       },
+      get_area = function(value) {
+        if (missing(value)) {
+          StochasticPatch___FF16w__FF16_Env__get_area__get(self)
+        } else {
+          stop("StochasticPatch<FF16w,FF16_Env>$get_area is read-only")
+        }
+      },
       parameters = function(value) {
         if (missing(value)) {
           StochasticPatch___FF16w__FF16_Env__parameters__get(self)
@@ -5171,6 +5227,13 @@ StochasticPatch <- function(T, E) {
           stop("StochasticPatch<FF16r,FF16_Env>$height_max is read-only")
         }
       },
+      get_area = function(value) {
+        if (missing(value)) {
+          StochasticPatch___FF16r__FF16_Env__get_area__get(self)
+        } else {
+          stop("StochasticPatch<FF16r,FF16_Env>$get_area is read-only")
+        }
+      },
       parameters = function(value) {
         if (missing(value)) {
           StochasticPatch___FF16r__FF16_Env__parameters__get(self)
@@ -5287,6 +5350,13 @@ StochasticPatch <- function(T, E) {
           stop("StochasticPatch<K93,K93_Env>$height_max is read-only")
         }
       },
+      get_area = function(value) {
+        if (missing(value)) {
+          StochasticPatch___K93__K93_Env__get_area__get(self)
+        } else {
+          stop("StochasticPatch<K93,K93_Env>$get_area is read-only")
+        }
+      },
       parameters = function(value) {
         if (missing(value)) {
           StochasticPatch___K93__K93_Env__parameters__get(self)

inst/RcppR6_classes.yml

@@ -456,6 +456,7 @@ Patch:
   active:
     time: {type: double, access: member}
     size: {type: size_t, access: member}
+    get_area: {type: double, access: member}
     height_max: {type: double, access: member}
     parameters: {type: "plant::Parameters<T,E>", access: member, name_cpp: r_parameters}
     environment: {type: E, access: member, name_cpp: r_environment}
@@ -652,6 +653,7 @@ StochasticPatch:
     time: {type: double, access: member}
     size: {type: size_t, access: member}
     height_max: {type: double, access: member}
+    get_area: {type: double, access: member}
     parameters: {type: "plant::Parameters<T,E>", access: member, name_cpp: r_parameters}
     environment: {type: E, access: member, name_cpp: r_environment}
     species: {type: "std::vector<plant::StochasticSpecies<T,E> >", access: member, name_cpp: r_species}

inst/include/plant/patch.h

@@ -28,7 +28,7 @@ class Patch {
 
   //Try using pointer in place of object itself
   double time() const {return environment.time;}
-
+  double get_area() const { return area;}
   double height_max() const;
 
   double compute_competition(double height) const;
@@ -131,6 +131,8 @@ class Patch {
   void compute_rates();
 
   parameters_type parameters;
+
+  double area;
   environment_type environment;
   std::vector<species_type> species;
 
@@ -145,6 +147,7 @@ class Patch {
 template <typename T, typename E>
 Patch<T,E>::Patch(parameters_type p, environment_type e, Control c)
   : parameters(p),
+    area(p.patch_area),
     environment(e),
     control(c),
     environment_cache(6) {  // length of ode::Step
@@ -221,7 +224,7 @@ double Patch<T,E>::compute_competition(double height) const {
   double tot = 0.0;
   for (size_t i = 0; i < species.size(); ++i) {
     if (!is_mutant_run) {
-      tot += species[i].compute_competition(height);
+      tot += species[i].compute_competition(height) / area;
     }
   }
   return tot;

inst/include/plant/scm.h

@@ -75,7 +75,7 @@ SCM<T, E>::SCM(parameters_type p, environment_type e, Control c)
 
   parameters.validate();
   if (!util::identical(parameters.patch_area, 1.0)) {
-    util::stop("Patch area must be exactly 1 for the SCM");
+    util::warning("We recommened keeping patch_area = 1 for the SCM, as need to check units for all other sizes");
   }
 }
 

inst/include/plant/stochastic_patch.h

@@ -26,6 +26,7 @@ class StochasticPatch {
 
   size_t size() const {return species.size();}
   double time() const {return environment.time;}
+  double get_area() const { return area;}
 
   double height_max() const;
 
@@ -44,6 +45,8 @@ class StochasticPatch {
   // * ODE interface
   size_t ode_size() const;
   double ode_time() const;
+  double area;
+
   ode::const_iterator set_ode_state(ode::const_iterator it, double time);
   ode::iterator       ode_state(ode::iterator it) const;
   ode::iterator       ode_rates(ode::iterator it) const;
@@ -85,7 +88,7 @@ class StochasticPatch {
 template <typename T, typename E>
 StochasticPatch<T,E>::StochasticPatch(parameters_type p, environment_type e, Control c)
   : parameters(p),
-
+    area(p.patch_area),
     environment(e),
     control(c) {
   parameters.validate();
@@ -119,7 +122,7 @@ template <typename T, typename E>
 double StochasticPatch<T,E>::compute_competition(double height) const {
   double tot = 0.0;
   for (size_t i = 0; i < species.size(); ++i) {
-      tot += species[i].compute_competition(height);
+    tot += species[i].compute_competition(height) / area;
   }
   return tot;
 }

inst/include/plant/util.h

@@ -112,6 +112,8 @@ bool is_decreasing(ForwardIterator first, ForwardIterator last) {
 
 void stop(const std::string&);
 
+void warning(const std::string &);
+
 template<typename T>
 std::string to_string(T x) {
   return std::to_string(x);