1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
|
% fallinthewind.mp
% L. Nobre G.
% 2014
input featpost3Dplus2D;
prologues := 1;
Spread := 9;
f := 4.5*(2.5,-1.3,1.7);
viewcentr := (0,0,2);
background := 0.9white;
color referencewind;
numeric dragcoef, dragcoeg, horangle, quadcoef;
horangle = 90;
referencewind = red*cosd( horangle ) + green*sind( horangle );
dragcoef = 0.065;
dragcoeg = 0.009;
quadcoef = 0.99;
def vecfunc( expr posit, veloc )=
begingroup
save heig, wind, field, wint, blow;
color wind, field, blow;
numeric heig, wint;
heig = Z( posit );
wind = referencewind*quadcoef*sqrt( heig+0.05 );
wind := wind-veloc;
wint = conorm( wind );
blow = N( wind );
field = (dragcoef+dragcoeg*wint)*wind-blue;
( field )
endgroup
enddef;
beginfig(1);
numeric i, j, k, numa, sa, ea, diffstep, ray;
color stp, svl, arrbe, arren;
path oneline[];
pen grossa, fina;
arrbe = (0,-8,5);
arren = (0,-5.5,5);
setthestage( 18, 20 );
TDAtiplen := 0.75;
TDAhalftipbase := 0.3;
TDAhalfthick := 0.15;
grossa = pencircle scaled 1.8pt;
fina = pencircle scaled 1.2pt;
diffstep = 0.015;
ray = 3.95;
numa = 3.45;
sa = 38;
ea = 142;
j = 0;
for i=sa step numa until ea:
j := j+1;
stp := black;
svl := ray*(green*cosd(i)+blue*sind(i));
oneline[j] := dragtrajectorypath( stp, svl, diffstep, vecfunc );
draw oneline[j] withpen grossa withcolor white;
endfor;
for k=1 upto j:
draw oneline[k] withpen fina withcolor red;
endfor;
tdarrow( arrbe, arren );
draw rp(black)--rp((0,0,4.5)) dashed withdots withpen fina withcolor blue;
produce_vga_border;
endfig;
end.
ray := 3.25;
numeric min, rst;
min = 0.5;
rst = 0.1;
for i=min step rst until ray:
%show i;
numeric xa, xb, xc, fa, fb, fc;
xa = 85;
svl := i*(green*cosd(xa)+blue*sind(xa));
oneline := dragtrajectorypath( black, svl, diffstep, vecfunc );
fa = PrintStep;
xb = 89;
svl := i*(green*cosd(xb)+blue*sind(xb));
oneline := dragtrajectorypath( black, svl, diffstep, vecfunc );
fb = PrintStep;
xc = 0;
forever:
exitif ( abs( fb - fa ) < 0.005 );
xc := xb - fb*(xb-xa)/(fb-fa);
svl := i*(green*cosd(xc)+blue*sind(xc));
oneline := dragtrajectorypath( black, svl, diffstep, vecfunc );
fc := PrintStep;
exitif ( abs( fc - fa ) < 0.0005 );
xa := xb;
fa := fb;
xb := xc;
fb := fc;
endfor;
message decimal( i ) & " " & decimal( xc );
endfor;
message " ";
|