Line data Source code
1 : /*
2 : This program is free software: you can redistribute it and/or modify
3 : it under the terms of the GNU General Public License as published by
4 : the Free Software Foundation, either version 3 of the License, or
5 : (at your option) any later version.
6 :
7 : This program is distributed in the hope that it will be useful,
8 : but WITHOUT ANY WARRANTY; without even the implied warranty of
9 : MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
10 : GNU General Public License for more details.
11 :
12 : You should have received a copy of the GNU General Public License
13 : along with this program. If not, see <http://www.gnu.org/licenses/>.
14 : */
15 : /*
16 : simple electric motor simulator class
17 : */
18 :
19 : #include "SIM_Motor.h"
20 : #include <AP_Motors/AP_Motors.h>
21 :
22 : using namespace SITL;
23 :
24 : // calculate rotational accel and thrust for a motor
25 452958054 : void Motor::calculate_forces(const struct sitl_input &input,
26 : uint8_t motor_offset,
27 : Vector3f &torque,
28 : Vector3f &thrust,
29 : const Vector3f &velocity_air_bf,
30 : const Vector3f &gyro,
31 : float air_density,
32 : float voltage,
33 : bool use_drag)
34 : {
35 :
36 452958054 : const float pwm = input.servos[motor_offset+servo];
37 : float command;
38 452958054 : float thrust_sign = 1;
39 452958054 : if (rsc_servo >= 0) {
40 : // variable-pitch rotor; the motor's servo channel commands
41 : // blade pitch about mid-PWM and rotor speed comes from the
42 : // shared rotor-speed-control servo. The governor is assumed
43 : // to hold the rotor at nominal speed from half RSC output;
44 : // thrust scales with blade pitch and the square of rotor
45 : // speed
46 2214400 : const float rotor_speed = constrain_float((input.servos[motor_offset+rsc_servo]-1000)*0.002, 0, 1);
47 2214400 : const float collective = constrain_float((pwm-1500)*0.002, -1, 1);
48 2214400 : if (is_negative(collective)) {
49 1244803 : thrust_sign = -1;
50 : }
51 2214400 : command = fabsf(collective) * sq(rotor_speed);
52 : } else {
53 450743654 : command = pwm_to_command(pwm);
54 : }
55 452958054 : float voltage_scale = voltage / voltage_max;
56 :
57 452958054 : if (voltage_scale < 0.1) {
58 : // battery is dead
59 4764 : torque.zero();
60 4764 : thrust.zero();
61 4764 : current = 0;
62 4764 : return;
63 : }
64 :
65 : // apply slew limiter to command
66 452953290 : uint64_t now_us = AP_HAL::micros64();
67 452953290 : if (last_calc_us != 0 && slew_max > 0) {
68 452419854 : float dt = (now_us - last_calc_us)*1.0e-6;
69 452419854 : float slew_max_change = slew_max * dt;
70 452419854 : command = constrain_float(command, last_command-slew_max_change, last_command+slew_max_change);
71 : }
72 452953290 : last_calc_us = now_us;
73 452953290 : last_command = command;
74 :
75 : // velocity of motor through air
76 452953290 : Vector3f motor_vel = velocity_air_bf;
77 :
78 : // add velocity of motor about center due to vehicle rotation
79 452953290 : motor_vel += -(position % gyro);
80 :
81 : // calculate velocity into prop, clipping at zero
82 452953290 : float velocity_in = MAX(0, -motor_vel.projected(thrust_vector).z);
83 :
84 : // get thrust for untilted motor
85 452953290 : float motor_thrust = thrust_sign * calc_thrust(command, air_density, velocity_in, voltage_scale);
86 :
87 : // the yaw torque of the motor. A variable-pitch rotor spins in
88 : // the same direction regardless of blade pitch sign, so the drag
89 : // torque reaction on the frame does not reverse with negative
90 : // thrust
91 452953290 : const float yaw_scale = 0.05 * diagonal_size * motor_thrust * thrust_sign;
92 452953290 : Vector3f rotor_torque = thrust_vector * yaw_factor * command * yaw_scale * -1.0;
93 :
94 : // thrust in bodyframe NED
95 452953290 : thrust = thrust_vector * motor_thrust;
96 :
97 : // work out roll and pitch of motor relative to it pointing straight up
98 452953290 : float roll = 0, pitch = 0;
99 :
100 452953290 : uint64_t now = AP_HAL::micros64();
101 :
102 : // possibly roll and/or pitch the motor
103 452953290 : if (roll_servo >= 0) {
104 1277918 : uint16_t servoval = update_servo(input.servos[roll_servo+motor_offset], now, last_roll_value);
105 1277918 : if (roll_min < roll_max) {
106 0 : roll = constrain_float(roll_min + (servoval-1000)*0.001*(roll_max-roll_min), roll_min, roll_max);
107 : } else {
108 1277918 : roll = constrain_float(roll_max + (2000-servoval)*0.001*(roll_min-roll_max), roll_max, roll_min);
109 : }
110 : }
111 452953290 : if (pitch_servo >= 0) {
112 6622906 : uint16_t servoval = update_servo(input.servos[pitch_servo+motor_offset], now, last_pitch_value);
113 6622906 : if (pitch_min < pitch_max) {
114 0 : pitch = constrain_float(pitch_min + (servoval-1000)*0.001*(pitch_max-pitch_min), pitch_min, pitch_max);
115 : } else {
116 6622906 : pitch = constrain_float(pitch_max + (2000-servoval)*0.001*(pitch_min-pitch_max), pitch_max, pitch_min);
117 : }
118 : }
119 452953290 : last_change_usec = now;
120 :
121 : // possibly rotate the thrust vector and the rotor torque
122 452953290 : if (!is_zero(roll) || !is_zero(pitch)) {
123 4845874 : Matrix3f rotation;
124 4845874 : rotation.from_euler(radians(roll), radians(pitch), 0);
125 4845874 : thrust = rotation * thrust;
126 4845874 : rotor_torque = rotation * rotor_torque;
127 : }
128 :
129 452953290 : if (use_drag) {
130 : // Momentum drag is modelled as isotropic, scaled by the total
131 : // thrust, which handles tilted motors correctly. The component
132 : // along the rotor axis is really a 'hack' to compensate for
133 : // incorrect modelling of the variation of thrust with inflow
134 : // velocity. If not applied, the vehicle will climb at an
135 : // unrealistic rate during operation in STABILIZE. TODO replace
136 : // prop and motor model with one based on DC motor, momentum
137 : // disc and blade element theory.
138 452953290 : const float momentum_drag_factor = momentum_drag_coefficient * sqrtf(air_density * true_prop_area * thrust.length());
139 452953290 : thrust -= motor_vel * momentum_drag_factor;
140 : }
141 :
142 : // calculate total torque in newton-meters
143 452953290 : torque = (position % thrust) + rotor_torque;
144 :
145 : // calculate current
146 452953290 : float power = power_factor * fabsf(motor_thrust);
147 452953290 : current = power / MAX(voltage, 0.1);
148 : }
149 :
150 : /*
151 : update and return current value of a servo. Calculated as 1000..2000
152 : */
153 7900824 : uint16_t Motor::update_servo(uint16_t demand, uint64_t time_usec, float &last_value) const
154 : {
155 7900824 : if (servo_rate <= 0) {
156 0 : return demand;
157 : }
158 7900824 : if (servo_type == SERVO_RETRACT) {
159 : // handle retract servos
160 731082 : if (demand > 1700) {
161 333956 : demand = 2000;
162 397126 : } else if (demand < 1300) {
163 397126 : demand = 1000;
164 : } else {
165 0 : demand = last_value;
166 : }
167 : }
168 7900824 : demand = constrain_int16(demand, 1000, 2000);
169 7900824 : float dt = (time_usec - last_change_usec) * 1.0e-6f;
170 : // assume servo moves through 90 degrees over 1000 to 2000
171 7900824 : float max_change = 1000 * (dt / servo_rate) * 60.0f / 90.0f;
172 7900824 : last_value = constrain_float(demand, last_value-max_change, last_value+max_change);
173 7900824 : return uint16_t(last_value+0.5);
174 : }
175 :
176 :
177 : // calculate current and voltage
178 452958054 : float Motor::get_current(void) const
179 : {
180 452958054 : return current;
181 : }
182 :
183 : // setup PWM ranges for this motor
184 452962818 : void Motor::setup_params(uint16_t _pwm_min, uint16_t _pwm_max, float _spin_min, float _spin_max, float _expo, float _slew_max,
185 : float _diagonal_size, float _power_factor, float _voltage_max, float _effective_prop_area,
186 : float _velocity_max, Vector3f _position, Vector3f _thrust_vector, float _yaw_factor,
187 : float _true_prop_area, float _momentum_drag_coefficient)
188 : {
189 452962818 : mot_pwm_min = _pwm_min;
190 452962818 : mot_pwm_max = _pwm_max;
191 452962818 : mot_spin_min = _spin_min;
192 452962818 : mot_spin_max = _spin_max;
193 452962818 : mot_expo = _expo;
194 452962818 : slew_max = _slew_max;
195 452962818 : power_factor = _power_factor;
196 452962818 : voltage_max = _voltage_max;
197 452962818 : effective_prop_area = _effective_prop_area;
198 452962818 : max_outflow_velocity = _velocity_max;
199 452962818 : true_prop_area = _true_prop_area;
200 452962818 : momentum_drag_coefficient = _momentum_drag_coefficient;
201 452962818 : diagonal_size = _diagonal_size;
202 :
203 452962818 : if (!_position.is_zero()) {
204 0 : position = _position;
205 : } else {
206 452962818 : position.x = cosf(radians(angle)) * _diagonal_size;
207 452962818 : position.y = sinf(radians(angle)) * _diagonal_size;
208 452962818 : position.z = 0;
209 : }
210 :
211 452962818 : if (!_thrust_vector.is_zero()) {
212 0 : thrust_vector = _thrust_vector;
213 : }
214 452962818 : if (!is_zero(_yaw_factor)) {
215 0 : yaw_factor = _yaw_factor;
216 : }
217 452962818 : }
218 :
219 : /*
220 : convert a PWM value to a command value from 0 to 1
221 : */
222 450743654 : float Motor::pwm_to_command(float pwm) const
223 : {
224 450743654 : const float pwm_thrust_max = mot_pwm_min + mot_spin_max * (mot_pwm_max - mot_pwm_min);
225 450743654 : const float pwm_thrust_min = mot_pwm_min + mot_spin_min * (mot_pwm_max - mot_pwm_min);
226 450743654 : const float pwm_thrust_range = pwm_thrust_max - pwm_thrust_min;
227 450743654 : return constrain_float((pwm-pwm_thrust_min)/pwm_thrust_range, 0, 1);
228 : }
229 :
230 : /*
231 : calculate thrust given a command value
232 : */
233 452953290 : float Motor::calc_thrust(float command, float air_density, float velocity_in, float voltage_scale) const
234 : {
235 452953290 : float velocity_out = voltage_scale * max_outflow_velocity * sqrtf((1-mot_expo)*command + mot_expo*sq(command));
236 452953290 : float ret = 0.5 * air_density * effective_prop_area * (sq(velocity_out) - sq(velocity_in));
237 : #if 0
238 : if (command > 0) {
239 : ::printf("air_density=%f effective_prop_area=%f velocity_in=%f velocity_max=%f\n",
240 : air_density, effective_prop_area, velocity_in, voltage_scale * max_outflow_velocity);
241 : ::printf("calc_thrust %.3f %.3f\n", command, ret);
242 : }
243 : #endif
244 452953290 : return ret;
245 : }
|