-
Notifications
You must be signed in to change notification settings - Fork 2
/
complex.fs
166 lines (131 loc) · 5.87 KB
/
complex.fs
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
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
\ complex numbers
\ Authors: Anton Ertl, Bernd Paysan
\ Copyright (C) 2005,2007,2015,2019,2020,2021,2022 Free Software Foundation, Inc.
\ This file is part of Gforth.
\ Gforth is free software; you can redistribute it and/or
\ modify it under the terms of the GNU General Public License
\ as published by the Free Software Foundation, either version 3
\ of the License, or (at your option) any later version.
\ This program is distributed in the hope that it will be useful,
\ but WITHOUT ANY WARRANTY; without even the implied warranty of
\ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
\ GNU General Public License for more details.
\ You should have received a copy of the GNU General Public License
\ along with this program. If not, see http://www.gnu.org/licenses/.
\ *** Complex arithmetic *** 23sep91py
: complex' ( n -- offset ) 2* floats ;
: complex+ ( zaddr -- zaddr' ) float+ float+ ;
\ simple operations 02mar05py
: fl> ( -- r ) f@local0 lp+ ;
: zdup ( z -- z z ) fover fover ;
: zdrop ( z -- ) fdrop fdrop ;
: zover ( z1 z2 -- z1 z2 z1 ) 3 fpick 3 fpick ;
: z>r ( z -- l:z ) f>l f>l ;
: zr> ( r:z -- z ) fl> fl> ;
: zswap ( z1 z2 -- z2 z1 ) frot f>l frot fl> ;
: zpick ( z1 .. zn n -- z1 .. zn z1 ) 2* 1+ dup >r fpick r> fpick ;
\ : zpin 2* 1+ dup >r fpin r> fpin ;
: zdepth ( -- u ) fdepth 2/ ;
: zrot ( z1 z2 z3 -- z2 z3 z1 ) z>r zswap zr> zswap ;
: z-rot ( z1 z2 z3 -- z3 z1 z2 ) zswap z>r zswap zr> ;
: z@ ( zaddr -- z ) dup >r f@ r> float+ f@ ;
: z! ( z zaddr -- ) dup >r float+ f! r> f! ;
: z+! ( z zaddr -- ) dup >r float+ f+! r> f+! ;
\ locals 10jan15py
Create z!-table ' z! , ' z+! ,
: to-z: ( -- ) -14 throw ;
to-opt: ( !!?addr!! ) POSTPONE laddr# >body @ lp-offset, z!-table to-!, ;
: compile-pushlocal-z ( a-addr -- ) ( run-time: z -- )
locals-size @ alignlp-f float+ float+ dup locals-size !
swap !
]] f>l f>l [[ ;
: compile-z@local ( n -- )
case
0 of ]] f@local0 f@local1 [[ endof
1 floats of ]] f@local1 f@local# [[ 2 floats , endof
dup postpone f@local# dup , postpone f@local# float+ ,
endcase ;
also locals-types definitions
: z: ( "name" -- a-addr xt )
create-local ['] to-z: set-to ['] compile-pushlocal-z
does> @ lp-offset compile-z@local ;
: z^ ( "name" -- a-addr xt )
w^ drop ['] compile-pushlocal-z ;
previous definitions
also locals-types
z: some-zlocal 2drop
previous
\ simple operations 02mar05py
: z+ ( z1 z2 -- z1+z2 ) frot f+ f>l f+ fl> ;
: z- ( z1 z2 -- z1-z2 ) fnegate frot f+ f>l f- fl> ;
: zr- ( z1 z2 -- z2-z1 ) frot f- f>l fswap f- fl> ;
: x+ ( z r -- z+r ) frot f+ fswap ;
: x- ( z r -- z-r ) fnegate x+ ;
: z* ( z1 z2 -- z1*z2 )
fdup 4 fpick f* f>l fover 3 fpick f* f>l
f>l fswap fl> f* f>l f* fl> f- fl> fl> f+ ;
: zscale ( z r -- z*r ) ftuck f* f>l f* fl> ;
\ simple operations 02mar05py
: znegate ( z -- -z ) fnegate fswap fnegate fswap ;
: zconj ( rr ri -- rr -ri ) fnegate ;
: z*i ( z -- z*i ) fnegate fswap ;
: z/i ( z -- z/i ) fswap fnegate ;
: zsqabs ( z -- |z|² ) fdup f* fswap fdup f* f+ ;
: 1/z ( z -- 1/z ) zconj zdup zsqabs 1/f zscale ;
: z/ ( z1 z2 -- z1/z2 ) 1/z z* ;
: |z| ( z -- r ) zsqabs fsqrt ;
: zabs ( z -- |z| ) |z| 0e ;
: z2/ ( z -- z/2 ) f2/ f>l f2/ fl> ;
: z2* ( z -- z*2 ) f2* f>l f2* fl> ;
: >polar ( z -- r theta ) zdup |z| fswap frot fatan2 ;
: polar> ( r theta -- z ) fsincos frot zscale fswap ;
\ zexp zln 02mar05py
: zexp ( z -- exp[z] ) fsincos fswap frot fexp zscale ;
: pln ( z -- pln[z] ) zdup fswap fatan2 frot frot |z| fln fswap ;
: zln ( z -- ln[z] ) >polar fswap fln fswap ;
: z0= ( z -- flag ) f0= >r f0= r> and ;
: zsqrt ( z -- sqrt[z] )
zdup z0= 0= IF
fdup f0= IF fdrop fsqrt 0e EXIT THEN
zln z2/ zexp THEN ;
: z** ( z1 z2 -- z1**z2 ) zswap zln z* zexp ;
\ Test: Fibonacci-Zahlen
1e 5e fsqrt f+ f2/ fconstant phi
: zfib ( z1 -- fib[z1] )
zdup z>r phi 0e zswap z**
zr> zswap z>r [ 1e phi f- ] FLiteral
0e zswap z** znegate zr> z+
[ 5e fsqrt 1/f ] FLiteral zscale ;
\ complexe Operationen 02mar05py
: zsinh ( z -- sinh[z] ) zexp zdup 1/z z- z2/ ;
: zcosh ( z -- cosh[z] ) zexp zdup 1/z z+ z2/ ;
: ztanh ( z -- tanh[z] ) z2* zexp zdup 1e 0e z- zswap 1e 0e z+ z/ ;
: zsin ( z -- sin[z] ) z*i zsinh z/i ;
: zcos ( z -- cos[z] ) z*i zcosh ;
: ztan ( z -- tan[z] ) z*i ztanh z/i ;
: Real ( z -- r ) fdrop ;
: Imag ( z -- i ) fnip ;
: Re ( z -- zr ) Real 0e ;
: Im ( z -- zi ) Imag 0e ;
\ complexe Operationen 02mar05py
: zasinh ( z -- asinh[z] ) zdup 1e f+ zover 1e f- z* zsqrt z+ pln ;
: zacosh ( z -- acosh[z] )
zdup 1e x- z2/ zsqrt zswap 1e x+ z2/ zsqrt z+ pln z2* ;
: zatanh ( z -- atanh[z] ) zdup 1e x+ zln zswap 1e x- znegate pln z- z2/ ;
: zacoth ( z -- acoth[z] ) znegate zdup 1e x- pln zswap 1e x+ pln z- z2/ ;
pi f2/ FConstant pi/2
: zasin ( z -- -iln[iz+sqrt[1-z^~2]] ) z*i zasinh z/i ;
: zacos ( z -- pi/2-asin[z] ) pi/2 0e zswap zasin z- ;
: zatan ( z -- [ln[1+iz]-ln[1-iz]]/2i ) z*i zatanh z/i ;
: zacot ( z -- [ln[[z+i]/[z-i]]/2i ) z*i zacoth z/i ;
\ Ausgabe 24sep05py
Defer fc. ' f. IS fc.
: z. ( z -- )
zdup z0= IF zdrop ." 0 " exit THEN
fdup f0= IF fdrop fc. exit THEN fswap
fdup f0= IF fdrop
ELSE fc. 1 backspaces
fdup f0> IF ." +" THEN THEN
fc. 1 backspaces ." i " ;
: z.s ( z1 .. zn -- z1 .. zn )
zdepth 0 ?DO i zpick zswap z>r z. zr> LOOP ;