rho_sw = 1025;
rho_air = 1.168;
ad = 0.47;  % drag area ... approximately equal to cable
cd = 1;  % drag coef 
ag = 1;  % generator area
effg = 0.26;  % generator efficiency (from a commercial model wind-> electrical)
ap = 0.2;  % propulsor area
effp = 0.5;  % propulsion efficiency

drag = 0.5 * rho_sw*ad*cd* c_speed.^2;

% From Matlab calcs...yuck...
vout = real(-(wind_speed)/3 + (4*(wind_speed).^2)./(3*(8*(wind_speed).^3 - 27*effg*(wind_speed).^3 + 3*sqrt(3)*(-16*effg*(wind_speed).^6 + 27*effg^2*(wind_speed).^6).^(1/2)).^(1/3)) + (8*(wind_speed).^3 - 27*effg*(wind_speed).^3 + 3*sqrt(3)*sqrt(-16*effg*(wind_speed).^6 + 27*effg.^2*(wind_speed).^6)).^(1/3)/3);
reactionforce = ag*rho_air*(wind_speed.^2-vout.^2)/2;

fu = drag .* u_20m ./ c_speed + reactionforce .* wind_u./wind_speed;
fv = drag .* v_20m ./ c_speed + reactionforce .* wind_v./wind_speed;
f = (fu.^2 + fv.^2).^0.5;

power = 0.5 * rho_sw*ap* (f/(rho_sw*ap)).^(3/2)/effp;