Drag as request option

src/client/js/otp/modules/planner/PlannerModule.js

@@ -364,6 +364,7 @@ otp.modules.planner.PlannerModule =
             if(this.showIntermediateStops) queryParams.showIntermediateStops = this.showIntermediateStops;
             if(this.watts) queryParams.watts = this.watts;
             if(this.weight) queryParams.weight = this.weight;
+            if(this.aerodynamicDrag) queryParams.aerodynamicDrag = this.aerodynamicDrag;
             if(this.minimumMicromobilitySpeed) queryParams.minimumMicromobilitySpeed = this.minimumMicromobilitySpeed;
             if(this.maximumMicromobilitySpeed) queryParams.maximumMicromobilitySpeed = this.maximumMicromobilitySpeed;
 

src/client/js/otp/widgets/tripoptions/TripOptionsWidget.js

@@ -1031,6 +1031,7 @@ otp.widgets.tripoptions.Micromobility =
 
         var html = '<div class="notDraggable">Watts: <input id="'+this.id+'-watts-value" type="text" style="width:30px;" value="250" />';
         html += '<div class="notDraggable">Weight (kg): <input id="'+this.id+'-weight-value" type="text" style="width:30px;" value="105" />';
+        html += '<div class="notDraggable">Aerodynamic Drag (Cd * A): <input id="'+this.id+'-drag-value" type="text" style="width:35px;" value="0.408" />';
         html += '<div class="notDraggable">Min Speed (m/s): <input id="'+this.id+'-minspeed-value" type="text" style="width:30px;" value="0.8" />';
         html += '<div class="notDraggable">Max Speed (m/s): <input id="'+this.id+'-maxspeed-value" type="text" style="width:30px;" value="12.5" />';
         html += "</div>"
@@ -1046,6 +1047,11 @@ otp.widgets.tripoptions.Micromobility =
                 watts : parseFloat($('#'+this_.id+'-watts-value').val()),
             });
         });
+        $('#'+this.id+'-drag-value').change(function() {
+            this_.tripWidget.inputChanged({
+                aerodynamicDrag : parseFloat($('#'+this_.id+'-drag-value').val()),
+            });
+        });
         $('#'+this.id+'-weight-value').change(function() {
             this_.tripWidget.inputChanged({
                 watts : parseFloat($('#'+this_.id+'-weight-value').val()),
@@ -1072,6 +1078,10 @@ otp.widgets.tripoptions.Micromobility =
         if(!isNaN(weightVal)) {
             $('#'+this.id+'-weight-value').val(weightVal);
         }
+        var dragVal = parseFloat(planData.queryParams.aerodynamicDrag)
+        if(!isNaN(dragVal)) {
+         $('#'+this.id+'-drag-value').val(dragVal);
+        }
         var minSpeedVal = parseFloat(planData.queryParams.minimumMicromobilitySpeed)
         if(!isNaN(minSpeedVal)) {
             $('#'+this.id+'-minSpeed-value').val(minSpeedVal);

src/main/java/org/opentripplanner/api/common/RoutingResource.java

@@ -449,6 +449,27 @@ public abstract class RoutingResource {
     @QueryParam("weight")
     private Double weight;
 
+    /**
+     * This is coefficient of drag and frontal area multiplied together. The equation for drag resistance and the
+     * extracted value is as follows:
+     *
+     * Fdrag = 0.5 * Cd * A * Rho * V^2
+     *               ⎣CdA_⎦
+     *
+     * See https://www.gribble.org/cycling/power_v_speed.html
+     *
+     * where
+     * Cd = coefficient of drag
+     * A = frontal area in m^2
+     * Rho = air density in kg / m^3
+     *
+     * You need a wind tunnel to properly measure the overall interaction of the coefficient of drag and frontal area,
+     * but a study showed that a comfortable bicycling position had a Cd * A value of 0.408.
+     * See https://www.cyclingpowerlab.com/CyclingAerodynamics.aspx
+     */
+    @QueryParam("aerodynamicDrag")
+    private Double aerodynamicDrag;
+
     /**
      * The minimum speed of a personal micromobility vehicle. This should only be used to avoid unreasonably slow times
      * on hills. If it is desired to model effectively impossible travel uphill (ie the vehicle can't reasonably be
@@ -818,6 +839,9 @@ protected RoutingRequest buildRequest() throws ParameterException {
         if (weight != null)
             request.weight = weight;
 
+        if (aerodynamicDrag != null)
+            request.aerodynamicDrag = aerodynamicDrag;
+
         if (minimumMicromobilitySpeed != null)
             request.minimumMicromobilitySpeed = minimumMicromobilitySpeed;
 

src/main/java/org/opentripplanner/routing/core/RoutingRequest.java

@@ -590,6 +590,26 @@ public class RoutingRequest implements Cloneable, Serializable {
      */
     public double weight = 105;
 
+    /**
+     * This is coefficient of drag and frontal area multiplied together. The equation for drag resistance and the
+     * extracted value is as follows:
+     *
+     * Fdrag = 0.5 * Cd * A * Rho * V^2
+     *               ⎣CdA_⎦
+     *
+     * See https://www.gribble.org/cycling/power_v_speed.html
+     *
+     * where
+     * Cd = coefficient of drag
+     * A = frontal area in m^2
+     * Rho = air density in kg / m^3
+     *
+     * You need a wind tunnel to properly measure the overall interaction of the coefficient of drag and frontal area,
+     * but a study showed that a comfortable bicycling position had a Cd * A value of 0.408.
+     * See https://www.cyclingpowerlab.com/CyclingAerodynamics.aspx
+     */
+    public double aerodynamicDrag = 0.408;
+
     /** Saves split edge which can be split on origin/destination search
      *
      * This is used so that TrivialPathException is thrown if origin and destination search would split the same edge

src/main/java/org/opentripplanner/routing/edgetype/StreetEdge.java

@@ -894,7 +894,8 @@ public double calculateSpeed(RoutingRequest options, TraverseMode traverseMode,
                     options.weight,
                     Math.atan(0), // 0 slope beta
                     getRollingResistanceCoefficient(),
-                    ElevationUtils.ZERO_ELEVATION_DRAG_RESISTIVE_FORCE_COMPONENT,
+                    options.aerodynamicDrag,
+                    ElevationUtils.ZERO_ELEVATION_AIR_DENSITY,
                     options.minimumMicromobilitySpeed,
                     options.maximumMicromobilitySpeed
                 ),
@@ -980,8 +981,7 @@ public double timeLowerBound(RoutingRequest options) {
      *              difficulty of traveling over bumpy roadways. This value is used to calculate the value of `Frg` as
      *              noted in the above equations.See this wikipedia page for a list of coefficients by various surface
      *              types: https://en.wikipedia.org/wiki/Rolling_resistance#Rolling_resistance_coefficient_examples
-     * @param aerodynamicDragComponent The product of the coefficient of aerodynamic drag, frontal area and air density.
-     *              This value is product of (Cd * A * ρ) as noted in the above mathematical equations.
+     * @param airDensity The air density .
      * @param minSpeed The minimum speed that the micromobility should travel at in cases where the slope is so steep
      *              that it would be faster to walk with the vehicle.
      * @param maxSpeed The maximum speed the vehicle can travel at.
@@ -992,7 +992,8 @@ public static double calculateMicromobilitySpeed(
         double weight,
         double beta,
         double coefficientOfRollingResistance,
-        double aerodynamicDragComponent,
+        double aerodynamicDrag,
+        double airDensity,
         double minSpeed,
         double maxSpeed
     ) {
@@ -1022,26 +1023,28 @@ public static double calculateMicromobilitySpeed(
             weight *
             // These cosine and sine calculations could be precalculated during graph build
             (coefficientOfRollingResistance * Math.cos(beta) + Math.sin(beta));
+        // The interaction of aerodynamics and air density (Cd * A * p)
+        double dragResistance = aerodynamicDrag * airDensity;
 
         double a = (
             -Math.pow(dynamicRollingResistance, 3) / 27.0
         ) + (
             (2.0 * normalizedRollingFriction * dynamicRollingResistance) /
-                (3.0 * Math.pow(aerodynamicDragComponent, 2))
+                (3.0 * Math.pow(dragResistance, 2))
         ) + (
-            watts / aerodynamicDragComponent
+            watts / (dragResistance)
         );
         double b = (
-            2.0 / (9.0 * aerodynamicDragComponent)
+            2.0 / (9.0 * dragResistance)
         ) * (
             3.0 * normalizedRollingFriction -
                 (
-                    (2.0 * dynamicRollingResistance) / aerodynamicDragComponent
+                    (2.0 * dynamicRollingResistance) / (dragResistance)
                 )
         );
 
         double cardanicCheck = Math.pow(a, 2) + Math.pow(b, 3);
-        double rollingDragComponent = 2.0 / 3.0 * dynamicRollingResistance / aerodynamicDragComponent;
+        double rollingDragComponent = 2.0 / 3.0 * dynamicRollingResistance / (dragResistance);
         double speed;
         if (cardanicCheck >= 0) {
             double cardanicCheckSqrt = Math.sqrt(cardanicCheck);

src/main/java/org/opentripplanner/routing/edgetype/StreetWithElevationEdge.java

@@ -37,12 +37,12 @@ public class StreetWithElevationEdge extends StreetEdge {
     // an array of the length in meters of the corresponding gradient at the same index
     private short[] gradientLengths;
 
-    // The maximum resistive drag force component along this StreetWithElevationEdge. The difference of this resistive
-    // drag force component is likely extremely small along the vast majority of edges in the graph. Therefore, don't
-    // store all values in an array like the gradients and gradient lengths. Instead, use the maximum resistive drag
-    // component which would correspond to drag resistive force at the minimum altitude seen on this edge. This is an
-    // overestimate of aerodynamic drag.
-    private double maximumDragResistiveForceComponent;
+    // The maximum air density seen along this StreetWithElevationEdge. The difference of this resistive drag force
+    // component is likely extremely small along the vast majority of edges in the graph. Therefore, don't store all
+    // values in an array like the gradients and gradient lengths. Instead, use the maximum air density which would
+    // correspond to the air density observed at the minimum altitude seen on this edge. This is an overestimate of air
+    // density.
+    private double maximumAirDensity;
 
     public StreetWithElevationEdge(StreetVertex v1, StreetVertex v2, LineString geometry,
             I18NString name, double length, StreetTraversalPermission permission, boolean back) {
@@ -80,7 +80,7 @@ public boolean setElevationProfile(PackedCoordinateSequence elev, boolean comput
 
         gradients = costs.gradients;
         gradientLengths = costs.gradientLengths;
-        maximumDragResistiveForceComponent = costs.maximumDragResistiveForceComponent;
+        maximumAirDensity = costs.maximumAirDensity;
 
         return costs.flattened;
     }
@@ -143,7 +143,8 @@ public double calculateSpeed(RoutingRequest options, TraverseMode traverseMode,
                     options.weight,
                     Math.atan(gradients[i] / 100.0),
                     this.getRollingResistanceCoefficient(),
-                    maximumDragResistiveForceComponent,
+                    options.aerodynamicDrag,
+                    maximumAirDensity,
                     options.minimumMicromobilitySpeed,
                     options.maximumMicromobilitySpeed
                 ),

src/main/java/org/opentripplanner/routing/util/ElevationUtils.java

@@ -43,28 +43,8 @@ public static double getDynamicRollingResistance(double beta) {
      */
     public static final double GRAVITATIONAL_ACCELERATION_CONSTANT = 9.80665;
 
-    // the Cd * A * p value at 0 elevation
-    public static final double ZERO_ELEVATION_DRAG_RESISTIVE_FORCE_COMPONENT = getDragResistiveForceComponent(0);
-
-    /**
-     * This is  coefficient of drag and frontal area multiplied together. The equation for drag resistance and the
-     * extracted value is as follows:
-     *
-     * Fdrag = 0.5 * Cd * A * Rho * V^2
-     *               ⎣CdA_⎦
-     *
-     * See https://www.gribble.org/cycling/power_v_speed.html
-     *
-     * where
-     * Cd = coefficient of drag
-     * A = frontal area in m^2
-     * Rho = air density in kg / m^3
-     *
-     * Apparently you need a wind tunnel to properly measure the coefficient of drag, so for now, assume the following:
-     * Cd = 0.63
-     * A = 0.6
-     */
-    private static final double FRONTAL_AREA_DRAG_COMPONENT = 0.63 * 0.6;
+    // the Rho value at 0 elevation
+    public static final double ZERO_ELEVATION_AIR_DENSITY = getAirDensity(0);
 
     // the air pressure at sea level in Pascals
     // see https://www.omnicalculator.com/physics/air-pressure-at-altitude
@@ -91,13 +71,10 @@ public static double getDynamicRollingResistance(double beta) {
     private static final double TEMPERATURE_DECLINE_PER_METER_OF_ELEVATION_GAIN = -9.8 / 1000;
 
     /**
-     * Calculates the components of drag resistance except for the velocity assuming travel through dry earthy air. The
-     * equation for drag resistance and the extracted value is as follows:
+     * Calculates the approximate air density for dry earthy air at a certain elevation. The value Rho is a part of
+     * the equation for drag resistance which is as follows:
      *
      * Fdrag = 0.5 * Cd * A * Rho * V^2
-     *             ⎣dragComponent⎦
-     *
-     * Note that the 0.5 is accounted for as a part of the final micromobility speed calcuations.
      *
      * See https://www.gribble.org/cycling/power_v_speed.html
      *
@@ -151,20 +128,18 @@ public static double getDynamicRollingResistance(double beta) {
      *
      * @param altitude The altitude in meters
      */
-    public static double getDragResistiveForceComponent(double altitude) {
+    public static double getAirDensity(double altitude) {
         double randomlyGuessedTemperature = A_RANDOM_OUTDOOR_TEMPERATURE_IN_KELVIN + (
             altitude * TEMPERATURE_DECLINE_PER_METER_OF_ELEVATION_GAIN
         );
-        return FRONTAL_AREA_DRAG_COMPONENT * (
-            AIR_PRESSURE_AT_SEA_LEVEL *
+        return AIR_PRESSURE_AT_SEA_LEVEL *
             Math.exp(
                 -GRAVITATIONAL_ACCELERATION_CONSTANT *
                     EARTHY_AIR_MOLAR_MASS *
                     altitude /
                     (UNIVERSAL_GAS_CONSTANT * randomlyGuessedTemperature)
             ) /
-            (SPECIFIC_GAS_CONSTANT_FOR_DRY_AIR * randomlyGuessedTemperature)
-        );
+            (SPECIFIC_GAS_CONSTANT_FOR_DRY_AIR * randomlyGuessedTemperature);
     }
 
     private static double[] getLengthsFromElevation(CoordinateSequence elev) {
@@ -213,7 +188,7 @@ public static SlopeCosts getSlopeCosts(CoordinateSequence elev, boolean slopeLim
                 false,
                 new byte[]{0},
                 new short[]{(short) trueLength},
-                getDragResistiveForceComponent(0)
+                getAirDensity(0)
             );
         }
         double lengthMultiplier = trueLength / flatLength;
@@ -236,7 +211,7 @@ public static SlopeCosts getSlopeCosts(CoordinateSequence elev, boolean slopeLim
 
             // add to existing gradient bin
             boolean gradientAdded = false;
-            double minCoordinatesAltitude = Math.min(coordinates[i + 1].x, coordinates[i].x);
+            double minCoordinatesAltitude = Math.min(coordinates[i + 1].y, coordinates[i].y);
             for (GradientBin bin : gradients) {
                 if (bin.gradient == iGradient) {
                     bin.distance += run;
@@ -288,14 +263,14 @@ public static SlopeCosts getSlopeCosts(CoordinateSequence elev, boolean slopeLim
         // convert gradient info into arrays of primitives
         byte[] gradientsArr = new byte[gradients.size()];
         short[] gradientLengthsArr = new short[gradients.size()];
-        double maximumDragResistiveForceComponent = Double.MIN_VALUE;
+        double maximumAirDensity = Double.MIN_VALUE;
         for (int i = 0; i < gradients.size(); i++) {
             GradientBin bin = gradients.get(i);
             gradientsArr[i] = (byte) bin.gradient;
             gradientLengthsArr[i] = (short) bin.distance;
-            double dragReistiveForceComponent = getDragResistiveForceComponent(bin.minAltitude);
-            if (dragReistiveForceComponent > maximumDragResistiveForceComponent) {
-                maximumDragResistiveForceComponent = dragReistiveForceComponent;
+            double airDensity = getAirDensity(bin.minAltitude);
+            if (airDensity > maximumAirDensity) {
+                maximumAirDensity = airDensity;
             }
         }
 
@@ -312,7 +287,7 @@ public static SlopeCosts getSlopeCosts(CoordinateSequence elev, boolean slopeLim
             flattened,
             gradientsArr,
             gradientLengthsArr,
-            maximumDragResistiveForceComponent
+            maximumAirDensity
         );
     }
 

src/main/java/org/opentripplanner/routing/util/SlopeCosts.java

@@ -9,11 +9,11 @@ public class SlopeCosts {
     public final double lengthMultiplier; // Multiplier to get true length based on flat (projected) length
     public final byte[] gradients; // array of gradients as percents
     public final short[] gradientLengths; // array of the length of each gradient in meters
-    public final double maximumDragResistiveForceComponent; // the maximum resistive drag force component along an edge
+    public final double maximumAirDensity; // the maximum resistive drag force component along an edge
 
     public SlopeCosts(double slopeSpeedFactor, double slopeWorkFactor, double slopeSafetyCost,
                       double maxSlope, double lengthMultiplier, boolean flattened, byte[] gradients,
-                      short[] gradientLengths, double maximumDragResistiveForceComponent) {
+                      short[] gradientLengths, double maximumAirDensity) {
         this.slopeSpeedFactor = slopeSpeedFactor;
         this.slopeWorkFactor = slopeWorkFactor;
         this.slopeSafetyCost = slopeSafetyCost;
@@ -22,6 +22,6 @@ public SlopeCosts(double slopeSpeedFactor, double slopeWorkFactor, double slopeS
         this.flattened = flattened;
         this.gradients = gradients;
         this.gradientLengths = gradientLengths;
-        this.maximumDragResistiveForceComponent = maximumDragResistiveForceComponent;
+        this.maximumAirDensity = maximumAirDensity;
     }
 }