Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Blame
Revision Log

m_real_polymorph.f90 @ 430

Last change on this file since 430 was 430, checked in by touze, 11 years ago
import madx-5.01.00
File size: 175.2 KB

Line
1	!The Full Polymorphic Package
2	!Copyright (C) Etienne Forest
3	! Taylor polymorphism at execution is based on an idea
4	! and C++ prototype developed by J. Bengtsson circa 1990
5
6	module polymorphic_taylor
7	use complex_taylor
8	implicit none
9	public
10	logical(lp),private,parameter::t=.true.,f=.false.
11	integer,private::NO,ND,ND2,NP,NDPT,NV !,lastmaster ! 2002.12.13
12	INTEGER,PRIVATE::NMAX_pol=100
13	! real(dp),PRIVATE::EPS=c_1d_6
14	! integer ent,exi
15	integer,private,parameter::m1=mmmmmm1,m2=mmmmmm2,m3=mmmmmm3,ms=mmmmmm4
16	integer,private,parameter:: m11=m1+msm1,m12=m1+msm2,m13=m1+ms*m3, &
17	m21=m2+msm1,m22=m2+msm2,m23=m2+ms*m3, &
18	m31=m3+msm1,m32=m3+msm2,m33=m3+ms*m3
19	logical(lp),private::old
20	private resetpoly,resetpolyn ,resetpoly0,resetpolyn0 ,allocpoly,allocpolyn,resetpoly_R,resetpoly_Rn
21	PRIVATE K_OPT,A_OPT,e_OPT,ke_OPT
22	private allocenv,allocenvn,resetenv,resetenvn
23	private equal,Dequaldacon,equaldacon,iequaldacon ,realEQUAL,singleequal
24	private taylorEQUAL,EQUALtaylor,complexreal_8
25	private real_8univ,univreal_8
26	PRIVATE real_8MAP,MAPREAL_8
27	private unaryADD,add,daddsc,dscadd,addsc,scadd,iaddsc,iscadd
28	private unarySUB,subs,dscsub,dsubsc,scsub,subsc,iscsub,isubsc
29	private mul,dmulsc,dscmul,mulsc,scmul,imulsc,iscmul
30	private div,ddivsc,dscdiv,divsc,scdiv,idivsc,iscdiv
31	private POW,POWR,POWR8
32	private dexpt,dcost,dcosdt,dsindt,dsint,dlogt,dsqrtt
33	private greaterthan,greatereq
34	private lessthan,dlessthansc,dsclessthan,lessthansc,sclessthan,ilessthansc,isclessthan
35	private igreatersc,iscgreater,dgreatersc,dscgreater,greatersc,scgreater
36	private dgreatereqsc,dscgreatereq,greatereqsc,scgreatereq,igreatereqsc,iscgreatereq
37	private abst,full_abst
38	private lesseq,dlesseqsc,dsclesseq,lesseqsc,sclesseq,ilesseqsc,isclesseq
39	private eq,deqsc,dsceq,eqsc,sceq,ieqsc,isceq
40	private neq,dneqsc,dscneq,neqsc,scneq,ineqsc,iscneq
41	PRIVATE GETCHARnd2,GETintnd2,getchar,GETint,GETORDER,CUTORDER
42	private real_8ENV,ENVreal_8 ,RADTAYLORENV_8, ENV_8RADTAYLOR ,beamENV_8,normal_p
43	private absoftdatan2dr,absoftdcosdr,absoftdsindr,absoftdtandr,absoftdatandr
44	!complex stuff
45	private datant,datanDt,datan2t,dasint,dacost,dtant,dtandt,datanht
46	private dcosht,dsinht,dtanht,SINX_XT,SINHX_XT,polymorpht
47	! PRIVATE EQUAL1D,EQUAL2D
48	! end complex stuff
49	private printpoly,printdouble,printsingle,dmulmapconcat
50	private line
51	character(120) line
52	!integer npol
53	!parameter (npol=20)
54	!INTEGER iasspol
55	!type(real_8) DUMMYpol,DUMpol(ndum)
56	!integer iasspol0
57	!private assp
58	PRIVATE SINH_HR
59	PRIVATE SIN_HR
60	! PROBE_8 STUFF
61	PRIVATE RADTAYLORprobe_8,beamprobe_8
62	real(dp) :: sinhx_x_min=1e-4_dp
63	real(dp) :: sinhx_x_minp=1.0_dp ! 1.e-9 !c_0_0001
64
65
66	INTERFACE assignment (=)
67	MODULE PROCEDURE EQUAL ! 2002.10.9
68	! MODULE PROCEDURE EQUAL1D ! 2004.7.10
69	! MODULE PROCEDURE EQUAL2D ! 2004.7.10
70	MODULE PROCEDURE complexreal_8
71	!zgoubi
72	MODULE PROCEDURE realEQUAL !
73	MODULE PROCEDURE singleequal
74	MODULE PROCEDURE taylorEQUAL !taylor=real_8
75	MODULE PROCEDURE EQUALtaylor !real_8= taylor ! here 2002.10.9
76	MODULE PROCEDURE Dequaldacon !
77	MODULE PROCEDURE equaldacon !
78	! For resetting
79	MODULE PROCEDURE Iequaldacon
80	MODULE PROCEDURE Iequaldaconn
81	! For universal_taylor
82	MODULE PROCEDURE real_8univ
83	MODULE PROCEDURE univreal_8 ! new sagan
84	! for damap
85	MODULE PROCEDURE MAPreal_8
86	MODULE PROCEDURE real_8MAP
87	! FOR RADIATION
88	MODULE PROCEDURE real_8ENV
89	MODULE PROCEDURE ENVreal_8
90	MODULE PROCEDURE RADTAYLORENV_8
91	MODULE PROCEDURE ENV_8RADTAYLOR
92
93	! allowing normalization of polymorphic arrays directly
94	MODULE PROCEDURE beamENV_8
95	MODULE PROCEDURE normal_p
96	MODULE PROCEDURE RADTAYLORprobe_8
97	MODULE PROCEDURE beamprobe_8
98	end INTERFACE
99
100
101	!@ <table border="4" cellspacing="1" bordercolor="#000000" id="AutoNumber2" width="400" height="135">
102	!@ <tr>
103	!@ <td width="77" height="20" align="center">
104	!@ <span style="text-transform: uppercase">
105	!@ <font face="Times New Roman" size="1">+</font></span></td>
106	!@ <td width="77" height="20" align="center">
107	!@ <span style="text-transform: uppercase">
108	!@ <font face="Times New Roman" size="1">REAL_8</font></span></td>
109	!@ <td width="77" height="20" align="center">
110	!@ <span style="text-transform: uppercase">
111	!@ <font face="Times New Roman" size="1">
112	!@ Real(dp)</font></span></td>
113	!@ <td width="78" height="20" align="center">
114	!@ <span style="text-transform: uppercase">
115	!@ <font face="Times New Roman" size="1">Real(sp)</font></span></td>
116	!@ <td width="56" height="20" align="center">
117	!@ <span style="text-transform: uppercase">
118	!@ <font face="Times New Roman" size="1">Integer</font></span></td>
119	!@ </tr>
120	!@ <tr>
121	!@ <td width="77" height="20" align="center">
122	!@ <span style="text-transform: uppercase">
123	!@ <font face="Times New Roman" size="1">REAL_8</font></span></td>
124	!@ <td width="77" height="20" align="center">
125	!@ <span style="text-transform: uppercase; font-weight:700">
126	!@ <font face="Times New Roman" size="1">
127	!@ <a style="text-decoration: none" href="m_real_polymorph.htm#ADD">add</a></font></span></td>
128	!@ <td width="77" height="20" align="center">
129	!@ <span style="text-transform: uppercase; font-weight:700">
130	!@ <font face="Times New Roman" size="1">
131	!@ <a style="text-decoration: none" href="m_real_polymorph.htm#DADDSC">daddsc</a></font></span></td>
132	!@ <td width="78" height="20" align="center"><b>
133	!@ <font size="1" face="Times New Roman">
134	!@ <a style="text-decoration: none" href="m_real_polymorph.htm#ADDSC">ADDSC</a></font></b></td>
135	!@ <td width="56" height="20" align="center"><b>
136	!@ <font size="1" face="Times New Roman">
137	!@ <a style="text-decoration: none" href="m_real_polymorph.htm#IADDSC">IADDSC</a></font></b></td>
138	!@ </tr>
139	!@ <tr>
140	!@ <td width="77" height="20" align="center">
141	!@ <span style="text-transform: uppercase">
142	!@ <font face="Times New Roman" size="1">
143	!@ Real(dp)</font></span></td>
144	!@ <td width="77" height="20" align="center">
145	!@ <span style="text-transform: uppercase; font-weight:700">
146	!@ <font face="Times New Roman" size="1">
147	!@ <a style="text-decoration: none" href="m_real_polymorph.htm#DSCADD">dscadd</a></font></span></td>
148	!@ <td width="77" height="20" align="center">
149	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
150	!@ <td width="78" height="20" align="center">
151	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
152	!@ <td width="56" height="20" align="center">
153	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
154	!@ </tr>
155	!@ <tr>
156	!@ <td width="77" height="20" align="center">
157	!@ <span style="text-transform: uppercase">
158	!@ <font face="Times New Roman" size="1">Real(sp)</font></span></td>
159	!@ <td width="77" height="20" align="center"><b>
160	!@ <font size="1" face="Times New Roman">
161	!@ <a style="text-decoration: none" href="m_real_polymorph#SCADD">SCADD</a></font></b></td>
162	!@ <td width="77" height="20" align="center">
163	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
164	!@ <td width="78" height="20" align="center">
165	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
166	!@ <td width="56" height="20" align="center">
167	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
168	!@ </tr>
169	!@ <tr>
170	!@ <td width="77" height="20" align="center">
171	!@ <span style="text-transform: uppercase">
172	!@ <font face="Times New Roman" size="1">Integer</font></span></td>
173	!@ <td width="77" height="20" align="center"><b>
174	!@ <font size="1" face="Times New Roman">
175	!@ <a style="text-decoration: none" href="m_real_polymorph.htm#ISCADD">ISCADD</a></font></b></td>
176	!@ <td width="77" height="20" align="center">
177	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
178	!@ <td width="78" height="20" align="center">
179	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
180	!@ <td width="56" height="20" align="center">
181	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
182	!@ </tr>
183	!@ </table>
184
185	INTERFACE OPERATOR (+)
186	MODULE PROCEDURE unaryADD !
187	MODULE PROCEDURE add !
188	MODULE PROCEDURE daddsc !
189	MODULE PROCEDURE dscadd !
190	MODULE PROCEDURE addsc !
191	MODULE PROCEDURE scadd !
192	MODULE PROCEDURE iaddsc !
193	MODULE PROCEDURE iscadd !
194	END INTERFACE
195
196	!@ <table border="4" cellspacing="1" bordercolor="#000000" id="AutoNumber1" width="400" height="135">
197	!@ <tr>
198	!@ <td width="77" height="20" align="center">
199	!@ <span style="text-transform: uppercase">
200	!@ <font face="Times New Roman" size="1">-</font></span></td>
201	!@ <td width="77" height="20" align="center">
202	!@ <span style="text-transform: uppercase">
203	!@ <font face="Times New Roman" size="1">REAL_8</font></span></td>
204	!@ <td width="77" height="20" align="center">
205	!@ <span style="text-transform: uppercase">
206	!@ <font face="Times New Roman" size="1">
207	!@ Real(dp)</font></span></td>
208	!@ <td width="78" height="20" align="center">
209	!@ <span style="text-transform: uppercase">
210	!@ <font face="Times New Roman" size="1">Real(sp)</font></span></td>
211	!@ <td width="56" height="20" align="center">
212	!@ <span style="text-transform: uppercase">
213	!@ <font face="Times New Roman" size="1">Integer</font></span></td>
214	!@ </tr>
215	!@ <tr>
216	!@ <td width="77" height="20" align="center">
217	!@ <span style="text-transform: uppercase">
218	!@ <font face="Times New Roman" size="1">REAL_8</font></span></td>
219	!@ <td width="77" height="20" align="center">
220	!@ <span style="text-transform: uppercase">
221	!@ <font face="Times New Roman" size="1">
222	!@ <a style="text-decoration: none; font-weight: 700" href="m_real_polymorph.htm#SUBS">SUBS</a></font></span></td>
223	!@ <td width="77" height="20" align="center">
224	!@ <span style="text-transform: uppercase">
225	!@ <font face="Times New Roman" size="1">
226	!@ <a style="text-decoration: none; font-weight: 700" href="m_real_polymorph.htm#DSUBSC">dSUBsc</a></font></span></td>
227	!@ <td width="78" height="20" align="center">
228	!@ <font size="1" face="Times New Roman">
229	!@ <a style="text-decoration: none; font-weight: 700" href="m_real_polymorph.htm#SUBSC">SUBSC</a></font></td>
230	!@ <td width="56" height="20" align="center">
231	!@ <font size="1" face="Times New Roman">
232	!@ <a style="text-decoration: none; font-weight: 700" href="m_real_polymorph.htm#ISUBSC">ISUBSC</a></font></td>
233	!@ </tr>
234	!@ <tr>
235	!@ <td width="77" height="20" align="center">
236	!@ <span style="text-transform: uppercase">
237	!@ <font face="Times New Roman" size="1">
238	!@ Real(dp)</font></span></td>
239	!@ <td width="77" height="20" align="center">
240	!@ <span style="text-transform: uppercase">
241	!@ <font face="Times New Roman" size="1">
242	!@ <a style="text-decoration: none; font-weight: 700" href="m_real_polymorph.htm#DSCSUB">dscSUB</a></font></span></td>
243	!@ <td width="77" height="20" align="center">
244	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
245	!@ <td width="78" height="20" align="center">
246	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
247	!@ <td width="56" height="20" align="center">
248	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
249	!@ </tr>
250	!@ <tr>
251	!@ <td width="77" height="20" align="center">
252	!@ <span style="text-transform: uppercase">
253	!@ <font face="Times New Roman" size="1">Real(sp)</font></span></td>
254	!@ <td width="77" height="20" align="center">
255	!@ <font size="1" face="Times New Roman">
256	!@ <a style="text-decoration: none; font-weight: 700" href="m_real_polymorph.htm#SCSUB">SCSUB</a></font></td>
257	!@ <td width="77" height="20" align="center">
258	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
259	!@ <td width="78" height="20" align="center">
260	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
261	!@ <td width="56" height="20" align="center">
262	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
263	!@ </tr>
264	!@ <tr>
265	!@ <td width="77" height="20" align="center">
266	!@ <span style="text-transform: uppercase">
267	!@ <font face="Times New Roman" size="1">Integer</font></span></td>
268	!@ <td width="77" height="20" align="center">
269	!@ <font size="1" face="Times New Roman">
270	!@ <a style="text-decoration: none; font-weight: 700" href="m_real_polymorph.htm#ISCSUB">ISCSUB</a></font></td>
271	!@ <td width="77" height="20" align="center">
272	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
273	!@ <td width="78" height="20" align="center">
274	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
275	!@ <td width="56" height="20" align="center">
276	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
277	!@ </tr>
278	!@ </table>
279
280	INTERFACE OPERATOR (-)
281	MODULE PROCEDURE unarySUB !
282	MODULE PROCEDURE subs !
283	MODULE PROCEDURE dscsub !
284	MODULE PROCEDURE dsubsc !
285	MODULE PROCEDURE subsc !
286	MODULE PROCEDURE scsub !
287	MODULE PROCEDURE isubsc !
288	MODULE PROCEDURE iscsub !
289	END INTERFACE
290
291	!@ <table border="4" cellspacing="1" bordercolor="#000000" id="AutoNumber1" width="400" height="134">
292	!@ <tr>
293	!@ <td width="77" height="20" align="center">
294	!@ <span style="text-transform: uppercase">
295	!@ <font face="Times New Roman" size="1">*</font></span></td>
296	!@ <td width="77" height="20" align="center">
297	!@ <span style="text-transform: uppercase">
298	!@ <font face="Times New Roman" size="1">REAL_8</font></span></td>
299	!@ <td width="77" height="20" align="center">
300	!@ <span style="text-transform: uppercase">
301	!@ <font face="Times New Roman" size="1">
302	!@ Real(dp)</font></span></td>
303	!@ <td width="78" height="20" align="center">
304	!@ <span style="text-transform: uppercase">
305	!@ <font face="Times New Roman" size="1">Real(sp)</font></span></td>
306	!@ <td width="56" height="20" align="center">
307	!@ <span style="text-transform: uppercase">
308	!@ <font face="Times New Roman" size="1">Integer</font></span></td>
309	!@ </tr>
310	!@ <tr>
311	!@ <td width="77" height="20" align="center">
312	!@ <span style="text-transform: uppercase">
313	!@ <font face="Times New Roman" size="1">REAL_8</font></span></td>
314	!@ <td width="77" height="20" align="center">
315	!@ <span style="text-transform: uppercase">
316	!@ <font face="Times New Roman" size="1">
317	!@ <a style="text-decoration: none; font-weight:700" href="m_real_polymorph.htm#MUL">MUL</a></font></span></td>
318	!@ <td width="77" height="20" align="center">
319	!@ <span style="text-transform: uppercase">
320	!@ <font face="Times New Roman" size="1">
321	!@ <a style="text-decoration: none; font-weight:700" href="m_real_polymorph.htm#DMULSC">dMULsc</a></font></span></td>
322	!@ <td width="78" height="20" align="center">
323	!@ <font size="1" face="Times New Roman">
324	!@ <a style="text-decoration: none; font-weight:700" href="m_real_polymorph.htm#MULSC">MULSC</a></font></td>
325	!@ <td width="56" height="20" align="center">
326	!@ <font size="1" face="Times New Roman">
327	!@ <a style="text-decoration: none; font-weight:700" href="m_real_polymorph.htm#IMULSC">IMULSC</a></font></td>
328	!@ </tr>
329	!@ <tr>
330	!@ <td width="77" height="19" align="center">
331	!@ <span style="text-transform: uppercase">
332	!@ <font face="Times New Roman" size="1">
333	!@ Real(dp)</font></span></td>
334	!@ <td width="77" height="19" align="center">
335	!@ <span style="text-transform: uppercase">
336	!@ <font face="Times New Roman" size="1">
337	!@ <a style="text-decoration: none; font-weight:700" href="m_real_polymorph.htm#DSCMUL">dscMUL</a></font></span></td>
338	!@ <td width="77" height="19" align="center">
339	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
340	!@ <td width="78" height="19" align="center">
341	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
342	!@ <td width="56" height="19" align="center">
343	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
344	!@ </tr>
345	!@ <tr>
346	!@ <td width="77" height="20" align="center">
347	!@ <span style="text-transform: uppercase">
348	!@ <font face="Times New Roman" size="1">Real(sp)</font></span></td>
349	!@ <td width="77" height="20" align="center">
350	!@ <font size="1" face="Times New Roman">
351	!@ <a style="text-decoration: none; font-weight:700" href="m_real_polymorph.htm#SCMUL">SCMUL</a></font></td>
352	!@ <td width="77" height="20" align="center">
353	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
354	!@ <td width="78" height="20" align="center">
355	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
356	!@ <td width="56" height="20" align="center">
357	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
358	!@ </tr>
359	!@ <tr>
360	!@ <td width="77" height="20" align="center">
361	!@ <span style="text-transform: uppercase">
362	!@ <font face="Times New Roman" size="1">Integer</font></span></td>
363	!@ <td width="77" height="20" align="center">
364	!@ <font size="1" face="Times New Roman">
365	!@ <a style="text-decoration: none; font-weight:700" href="m_real_polymorph.htm#ISCMUL">ISCMUL</a></font></td>
366	!@ <td width="77" height="20" align="center">
367	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
368	!@ <td width="78" height="20" align="center">
369	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
370	!@ <td width="56" height="20" align="center">
371	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
372	!@ </tr>
373	!@ </table>
374
375	INTERFACE OPERATOR (*)
376	MODULE PROCEDURE mul !
377	MODULE PROCEDURE dmulsc !
378	MODULE PROCEDURE dscmul !
379	MODULE PROCEDURE mulsc !
380	MODULE PROCEDURE scmul !
381	MODULE PROCEDURE imulsc !
382	MODULE PROCEDURE iscmul !
383	MODULE PROCEDURE dmulmapconcat ! concatenation real_8*damap
384	END INTERFACE
385
386	!@ <table border="4" cellspacing="1" bordercolor="#000000" id="AutoNumber1" width="400" height="135">
387	!@ <tr>
388	!@ <td width="0" height="20" align="center">
389	!@ <span style="text-transform: uppercase">
390	!@ <font face="Times New Roman" size="1">/</font></span></td>
391	!@ <td width="0" height="20" align="center">
392	!@ <span style="text-transform: uppercase">
393	!@ <font face="Times New Roman" size="1">REAL_8</font></span></td>
394	!@ <td width="0" height="20" align="center">
395	!@ <span style="text-transform: uppercase">
396	!@ <font face="Times New Roman" size="1">
397	!@ Real(dp)</font></span></td>
398	!@ <td width="0" height="20" align="center">
399	!@ <span style="text-transform: uppercase">
400	!@ <font face="Times New Roman" size="1">Real(sp)</font></span></td>
401	!@ <td width="0" height="20" align="center">
402	!@ <span style="text-transform: uppercase">
403	!@ <font face="Times New Roman" size="1">Integer</font></span></td>
404	!@ </tr>
405	!@ <tr>
406	!@ <td width="0" height="20" align="center">
407	!@ <span style="text-transform: uppercase">
408	!@ <font face="Times New Roman" size="1">REAL_8</font></span></td>
409	!@ <td width="0" height="20" align="center">
410	!@ <span style="text-transform: uppercase">
411	!@ <font face="Times New Roman" size="1">
412	!@ <a style="text-decoration: none; font-weight: 700" href="m_real_polymorph.htm#DIV">div</a></font></span></td>
413	!@ <td width="0" height="20" align="center">
414	!@ <span style="text-transform: uppercase">
415	!@ <font face="Times New Roman" size="1">
416	!@ <a style="text-decoration: none; font-weight: 700" href="m_real_polymorph.htm#DDIVSC">dDIVsc</a></font></span></td>
417	!@ <td width="0" height="20" align="center"><font size="1">
418	!@ <a style="text-decoration: none; font-weight: 700" href="m_real_polymorph.htm#DIVSC">DIVSC</a></font></td>
419	!@ <td width="0" height="20" align="center"><font size="1">
420	!@ <a style="text-decoration: none; font-weight: 700" href="m_real_polymorph.htm#IDIVSC">IDIVSC</a></font></td>
421	!@ </tr>
422	!@ <tr>
423	!@ <td width="0" height="20" align="center">
424	!@ <span style="text-transform: uppercase">
425	!@ <font face="Times New Roman" size="1">
426	!@ Real(dp)</font></span></td>
427	!@ <td width="0" height="20" align="center">
428	!@ <span style="text-transform: uppercase">
429	!@ <font face="Times New Roman" size="1">
430	!@ <a style="text-decoration: none; font-weight: 700" href="m_real_polymorph.htm#DSCDIV">dscDIV</a></font></span></td>
431	!@ <td width="0" height="20" align="center">
432	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
433	!@ <td width="0" height="20" align="center">
434	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
435	!@ <td width="0" height="20" align="center">
436	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
437	!@ </tr>
438	!@ <tr>
439	!@ <td width="0" height="20" align="center">
440	!@ <span style="text-transform: uppercase">
441	!@ <font face="Times New Roman" size="1">Real(sp)</font></span></td>
442	!@ <td width="0" height="20" align="center"><font size="1">
443	!@ <a style="text-decoration: none; font-weight: 700" href="m_real_polymorph.htm#SCDIV">SCDIV</a></font></td>
444	!@ <td width="0" height="20" align="center">
445	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
446	!@ <td width="0" height="20" align="center">
447	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
448	!@ <td width="0" height="20" align="center">
449	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
450	!@ </tr>
451	!@ <tr>
452	!@ <td width="0" height="20" align="center">
453	!@ <span style="text-transform: uppercase">
454	!@ <font face="Times New Roman" size="1">Integer</font></span></td>
455	!@ <td width="0" height="20" align="center"><font size="1">
456	!@ <a style="text-decoration: none; font-weight: 700" href="m_real_polymorph.htm#ISCDIV">ISCDIV</a></font></td>
457	!@ <td width="0" height="20" align="center">
458	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
459	!@ <td width="0" height="20" align="center">
460	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
461	!@ <td width="0" height="20" align="center">
462	!@ <font size="2" face="Times New Roman"><b>F90</b></font></td>
463	!@ </tr>
464	!@ </table>
465
466
467
468	INTERFACE OPERATOR (/)
469	MODULE PROCEDURE div !
470	MODULE PROCEDURE ddivsc !
471	MODULE PROCEDURE dscdiv !
472	MODULE PROCEDURE divsc !
473	MODULE PROCEDURE scdiv !
474	MODULE PROCEDURE idivsc !
475	MODULE PROCEDURE iscdiv !
476	END INTERFACE
477
478	INTERFACE OPERATOR (**)
479	MODULE PROCEDURE POW !
480	MODULE PROCEDURE POWR !
481	MODULE PROCEDURE POWR8 !
482	END INTERFACE
483
484	! Logical Operators
485
486	INTERFACE OPERATOR (<)
487	MODULE PROCEDURE lessthan
488	MODULE PROCEDURE dlessthansc
489	MODULE PROCEDURE dsclessthan
490	MODULE PROCEDURE lessthansc
491	MODULE PROCEDURE sclessthan
492	MODULE PROCEDURE ilessthansc
493	MODULE PROCEDURE isclessthan
494	END INTERFACE
495
496	INTERFACE OPERATOR (>)
497	MODULE PROCEDURE greaterthan
498	MODULE PROCEDURE dgreatersc
499	MODULE PROCEDURE dscgreater
500	MODULE PROCEDURE greatersc
501	MODULE PROCEDURE scgreater
502	MODULE PROCEDURE igreatersc
503	MODULE PROCEDURE iscgreater
504	END INTERFACE
505
506	INTERFACE OPERATOR (>=)
507	MODULE PROCEDURE greatereq
508	MODULE PROCEDURE dgreatereqsc
509	MODULE PROCEDURE dscgreatereq
510	MODULE PROCEDURE greatereqsc
511	MODULE PROCEDURE scgreatereq
512	MODULE PROCEDURE igreatereqsc
513	MODULE PROCEDURE iscgreatereq
514	END INTERFACE
515
516	INTERFACE OPERATOR (<=)
517	MODULE PROCEDURE lesseq
518	MODULE PROCEDURE dlesseqsc
519	MODULE PROCEDURE dsclesseq
520	MODULE PROCEDURE lesseqsc
521	MODULE PROCEDURE sclesseq
522	MODULE PROCEDURE ilesseqsc
523	MODULE PROCEDURE isclesseq
524	END INTERFACE
525
526	INTERFACE OPERATOR (==)
527	MODULE PROCEDURE eq
528	MODULE PROCEDURE deqsc
529	MODULE PROCEDURE dsceq
530	MODULE PROCEDURE eqsc
531	MODULE PROCEDURE sceq
532	MODULE PROCEDURE ieqsc
533	MODULE PROCEDURE isceq
534	END INTERFACE
535
536	INTERFACE OPERATOR (/=)
537	MODULE PROCEDURE neq
538	MODULE PROCEDURE dneqsc
539	MODULE PROCEDURE dscneq
540	MODULE PROCEDURE neqsc
541	MODULE PROCEDURE scneq
542	MODULE PROCEDURE ineqsc
543	MODULE PROCEDURE iscneq
544	END INTERFACE
545
546
547
548	! New Operators
549
550	INTERFACE OPERATOR (.SUB.)
551	MODULE PROCEDURE GETORDER
552	MODULE PROCEDURE getchar
553	MODULE PROCEDURE GETint
554	END INTERFACE
555
556	INTERFACE OPERATOR (.CUT.)
557	MODULE PROCEDURE CUTORDER
558	END INTERFACE
559
560	INTERFACE OPERATOR (.PAR.)
561	MODULE PROCEDURE getcharnd2
562	MODULE PROCEDURE GETintnd2
563	END INTERFACE
564
565	! Intrinsic function (some are Fortran extension not always supported)
566
567	INTERFACE exp
568	MODULE PROCEDURE dexpt
569	END INTERFACE
570	INTERFACE dexp
571	MODULE PROCEDURE dexpt
572	END INTERFACE
573	INTERFACE cexp
574	MODULE PROCEDURE dexpt
575	END INTERFACE
576	INTERFACE cdexp
577	MODULE PROCEDURE dexpt
578	END INTERFACE
579
580	INTERFACE cos
581	MODULE PROCEDURE dcost
582	END INTERFACE
583	INTERFACE ccos
584	MODULE PROCEDURE dcost
585	END INTERFACE
586	INTERFACE cdcos
587	MODULE PROCEDURE dcost
588	END INTERFACE
589	INTERFACE dcos
590	MODULE PROCEDURE dcost
591	END INTERFACE
592
593	INTERFACE atanh
594	MODULE PROCEDURE datanht
595	END INTERFACE
596
597	INTERFACE tan
598	MODULE PROCEDURE dtant
599	END INTERFACE
600
601
602	INTERFACE dtan
603	MODULE PROCEDURE dtant
604	END INTERFACE
605
606	INTERFACE tand
607	MODULE PROCEDURE dtandt
608	! MODULE PROCEDURE absoftdtandr ! for non PC
609	END INTERFACE
610	INTERFACE dtand
611	MODULE PROCEDURE dtandt
612	! MODULE PROCEDURE absoftdtandr ! for non PC absoft
613	END INTERFACE
614
615
616	INTERFACE cosd
617	MODULE PROCEDURE dcosdt
618	! MODULE PROCEDURE absoftdcosdr ! for non PC
619	END INTERFACE
620	INTERFACE dcosd
621	MODULE PROCEDURE dcosdt
622	! MODULE PROCEDURE absoftdcosdr ! for non PC absoft
623	END INTERFACE
624
625
626	INTERFACE sin
627	MODULE PROCEDURE dsint
628	END INTERFACE
629	INTERFACE cdsin
630	MODULE PROCEDURE dsint
631	END INTERFACE
632	INTERFACE ccsin
633	MODULE PROCEDURE dsint
634	END INTERFACE
635	INTERFACE dsin
636	MODULE PROCEDURE dsint
637	END INTERFACE
638
639	INTERFACE sind
640	MODULE PROCEDURE dsindt
641	! MODULE PROCEDURE absoftdsindr ! for non PC
642	END INTERFACE
643	INTERFACE dsind
644	MODULE PROCEDURE dsindt
645	! MODULE PROCEDURE absoftdsindr ! for non PC absoft
646	END INTERFACE
647
648
649	INTERFACE log
650	MODULE PROCEDURE dlogt
651	END INTERFACE
652	INTERFACE dlog
653	MODULE PROCEDURE dlogt
654	END INTERFACE
655	INTERFACE cdlog
656	MODULE PROCEDURE dlogt
657	END INTERFACE
658	INTERFACE clog
659	MODULE PROCEDURE dlogt
660	END INTERFACE
661
662	INTERFACE sqrt
663	MODULE PROCEDURE dsqrtt
664	END INTERFACE
665	INTERFACE dsqrt
666	MODULE PROCEDURE dsqrtt
667	END INTERFACE
668
669	INTERFACE abs
670	MODULE PROCEDURE abst
671	END INTERFACE
672	INTERFACE dabs
673	MODULE PROCEDURE abst
674	END INTERFACE
675
676	! complex stuff
677	INTERFACE atan
678	MODULE PROCEDURE datant
679	END INTERFACE
680	INTERFACE datan
681	MODULE PROCEDURE datant
682	END INTERFACE
683
684	INTERFACE atan2
685	MODULE PROCEDURE datan2t
686	END INTERFACE
687	INTERFACE datan2
688	MODULE PROCEDURE datan2t
689	END INTERFACE
690
691	INTERFACE atan2d
692	MODULE PROCEDURE datan2dt
693	! MODULE PROCEDURE absoftdatan2dr ! for non PC
694	END INTERFACE
695	INTERFACE datan2d
696	MODULE PROCEDURE datan2dt
697	! MODULE PROCEDURE absoftdatan2dr ! for non PC absoft
698	END INTERFACE
699
700	INTERFACE atand
701	MODULE PROCEDURE datandt
702	! MODULE PROCEDURE absoftdatandr ! for non PC
703	END INTERFACE
704	INTERFACE datand
705	MODULE PROCEDURE datandt
706	! MODULE PROCEDURE absoftdatandr ! for non PC absoft
707	END INTERFACE
708
709	INTERFACE asin
710	MODULE PROCEDURE dasint
711	END INTERFACE
712	INTERFACE dasin
713	MODULE PROCEDURE dasint
714	END INTERFACE
715
716	INTERFACE acos
717	MODULE PROCEDURE dacost
718	END INTERFACE
719	INTERFACE dacos
720	MODULE PROCEDURE dacost
721	END INTERFACE
722
723	INTERFACE cosh
724	MODULE PROCEDURE dcosht
725	END INTERFACE
726	INTERFACE dcosh
727	MODULE PROCEDURE dcosht
728	END INTERFACE
729
730	INTERFACE sinh
731	MODULE PROCEDURE dsinht
732	END INTERFACE
733	INTERFACE dsinh
734	MODULE PROCEDURE dsinht
735	END INTERFACE
736
737	INTERFACE tanh
738	MODULE PROCEDURE dtanht
739	END INTERFACE
740	INTERFACE dtanh
741	MODULE PROCEDURE dtanht
742	END INTERFACE
743
744	! end Intrinsic function
745
746	! Non Intrinsic function
747
748	INTERFACE full_abs
749	MODULE PROCEDURE full_abst
750	END INTERFACE
751
752	INTERFACE morph
753	MODULE PROCEDURE polymorpht
754	END INTERFACE
755
756	INTERFACE SINHX_X
757	MODULE PROCEDURE SINH_HR
758	MODULE PROCEDURE SINHX_XT
759	END INTERFACE
760
761	INTERFACE SINX_X
762	MODULE PROCEDURE SIN_HR
763	MODULE PROCEDURE SINX_XT
764	END INTERFACE
765
766	! i/o
767
768	INTERFACE daprint
769	MODULE PROCEDURE printpoly !
770	MODULE PROCEDURE printdouble
771	MODULE PROCEDURE printsingle
772	END INTERFACE
773	INTERFACE print
774	MODULE PROCEDURE printpoly !
775	MODULE PROCEDURE printdouble
776	MODULE PROCEDURE printsingle
777	END INTERFACE
778
779
780
781	! constructors and destructors
782
783	INTERFACE alloc
784	MODULE PROCEDURE A_OPT !allocpoly !
785	MODULE PROCEDURE allocpolyn !
786	!radiation
787	MODULE PROCEDURE e_opt !
788	! MODULE PROCEDURE allocenv !
789	MODULE PROCEDURE allocenvn !
790	END INTERFACE
791
792	INTERFACE kill
793	MODULE PROCEDURE K_OPT !resetpoly0 !
794	MODULE PROCEDURE Ke_OPT !resetpoly0 !
795	MODULE PROCEDURE resetpolyn0 !
796	MODULE PROCEDURE resetpoly_R
797	MODULE PROCEDURE resetpoly_RN
798	!radiation
799	! MODULE PROCEDURE resetenv !
800	MODULE PROCEDURE resetenvn !
801	END INTERFACE
802
803	INTERFACE reset
804	MODULE PROCEDURE resetpoly !
805	MODULE PROCEDURE resetpolyn !
806	!radiation
807	MODULE PROCEDURE resetenv !
808	MODULE PROCEDURE resetenvn !
809	END INTERFACE
810
811	! Managing
812
813	INTERFACE ass
814	MODULE PROCEDURE assp
815	END INTERFACE
816
817	contains
818
819	subroutine make_it_knob(k,i,s)
820	implicit none
821	TYPE (real_8), intent(inout) :: k
822	real(dp), optional :: s
823	integer, intent(in) :: i
824	k%s=1.0_dp
825	if(present(s)) k%s=s
826	k%i=i
827	k%kind=3
828	end subroutine make_it_knob
829
830	subroutine kill_knob(k)
831	implicit none
832	TYPE (real_8), intent(inout) :: k
833	k%s=1.0_dp
834	k%i=0
835	k%kind=1
836	end subroutine kill_knob
837
838
839	FUNCTION polymorpht( S1 )
840	implicit none
841	TYPE (real_8) polymorpht
842	TYPE (taylor), INTENT (IN) :: S1
843	integer localmaster
844
845	localmaster=master
846	call ass(polymorpht)
847	polymorpht%t= s1
848	master=localmaster
849
850	END FUNCTION polymorpht
851
852
853	FUNCTION GETCHARnd2( S1, S2 )
854	implicit none
855	TYPE (real_8) GETCHARnd2
856	TYPE (real_8), INTENT (IN) :: S1
857	CHARACTER(*) , INTENT (IN) :: S2
858	integer localmaster
859	type(taylor) t
860
861	localmaster=master
862	call ass(GETCHARnd2)
863	call alloc(t)
864	t=s1
865
866	t=t.par.s2
867	GETCHARnd2%t=t
868
869
870
871	! GETCHARnd2%kind=m2
872	! if(s1%kind==m2) then
873	! GETCHARnd2%t=s1%t.par.s2
874	! else
875	! GETCHARnd2%t=zero
876	! endif
877
878	call kill(t)
879	master=localmaster
880
881	END FUNCTION GETCHARnd2
882
883	FUNCTION GETintnd2( S1, S2 )
884	implicit none
885	TYPE (real_8) GETintnd2
886	TYPE (real_8), INTENT (IN) :: S1
887	integer, INTENT (IN) :: S2(:)
888	integer localmaster
889	type(taylor) t
890
891	localmaster=master
892	call ass(GETintnd2)
893	call alloc(t)
894	t=s1
895	! if(s1%kind==m2) then
896	! GETintnd2%t=s1%t.par.s2
897	! else
898	! GETintnd2%kind=m2
899	! GETintnd2%t=zero
900	!
901	! endif
902	t=t.par.s2
903	GETintnd2%t=t
904
905
906	call kill(t)
907
908	master=localmaster
909
910	END FUNCTION GETintnd2
911
912	FUNCTION GETchar( S1, S2 )
913	implicit none
914	real(dp) GETchar
915	TYPE (real_8), INTENT (IN) :: S1
916	CHARACTER(*) , INTENT (IN) :: S2
917	! integer localmaster
918
919	GETchar=0.0_dp
920	if(s1%kind==m2) then
921	! GETchar%t=s1%t.sub.s2 ! OLD
922	GETchar=s1%t.sub.s2 ! CHANGE
923
924	endif
925
926	END FUNCTION GETchar
927
928	FUNCTION GETint( S1, S2 )
929	implicit none
930	real(dp) GETint
931	TYPE (real_8), INTENT (IN) :: S1
932	integer , INTENT (IN) :: S2(:)
933	integer i
934
935	GETint=0.0_dp
936	if(s1%kind==m2) then
937	GETint=s1%t.sub.s2 ! CHANGE
938	elseif(s1%kind==m1) then
939	!zgoubi
940	GETint=s1
941	do i=1,size(s2)
942	if(S2(i)/=0) then
943	GETint=0.0_dp
944	exit
945	endif
946	enddo
947
948	endif
949
950	END FUNCTION GETint
951
952	FUNCTION GETORDER( S1, S2 )
953	implicit none
954	TYPE (real_8) GETORDER
955	TYPE (real_8), INTENT (IN) :: S1
956	integer , INTENT (IN) :: S2
957	integer localmaster
958	localmaster=master
959	call ass(GETORDER)
960
961	if(s1%kind==m2) then
962	GETORDER%t=s1%t.sub.s2
963	else
964	GETORDER%kind=m1
965	GETORDER%r=0.0_dp
966	if(s2==0) GETORDER%r=s1%r
967	endif
968
969	master=localmaster
970
971	END FUNCTION GETORDER
972
973
974
975
976	FUNCTION CUTORDER( S1, S2 )
977	implicit none
978	TYPE (real_8) CUTORDER
979	TYPE (real_8), INTENT (IN) :: S1
980	integer , INTENT (IN) :: S2
981	integer localmaster
982
983	localmaster=master
984	call ass(CUTORDER)
985
986	CUTORDER=0.0_dp
987	if(s1%kind==m2) then
988	cutorder%kind=m2
989	CUTORDER%t=s1%t.CUT.s2
990	elseif(s1%kind==m1) then
991	if(s2>=1) CUTORDER%r=s1%r
992	cutorder%kind=m1
993	endif
994
995	master=localmaster
996
997	END FUNCTION CUTORDER
998
999
1000
1001
1002	FUNCTION greatereq( S1, S2 )
1003	implicit none
1004	logical(lp) greatereq
1005	TYPE (real_8), INTENT (IN) :: S1,S2
1006
1007	greatereq=.FALSE.
1008	select case(s1%kind+ms*s2%kind)
1009	case(m11)
1010	IF(S1%r>=S2%r) greatereq=.TRUE.
1011	case(m12,m21,m22)
1012	select case(s1%kind+ms*s2%kind)
1013	case(m21)
1014	IF((S1%t.SUB.'0')>=S2%r) greatereq=.TRUE.
1015	case(m12)
1016	IF(S1%r>=(S2%t.SUB.'0')) greatereq=.TRUE.
1017	case(m22)
1018	IF((S1%t.SUB.'0')>=(S2%t.SUB.'0')) greatereq=.TRUE.
1019	end select
1020	case(m13,m31,m32,m23,m33)
1021	select case(s1%kind+ms*s2%kind)
1022	case(m31,m13,m33)
1023	IF(S1%r>=S2%r) greatereq=.TRUE.
1024	case(m32)
1025	IF(S1%r>=(S2%t.SUB.'0')) greatereq=.TRUE.
1026	case(m23)
1027	IF((S1%t.SUB.'0')>=S2%r) greatereq=.TRUE.
1028	end select
1029	case default
1030	w_p=0
1031	w_p%nc=2
1032	w_p%fc='((1X,A72,/,1x,a72))'
1033	w_p%fi='(2((1X,i4)))'
1034	w_p%c(1)= " trouble in greatereq "
1035	w_p%c(2)= "s1%kind ,s2%kind "
1036	w_p=(/s1%kind ,s2%kind/)
1037	! call !write_e(0)
1038	end select
1039
1040	END FUNCTION greatereq
1041
1042	FUNCTION dgreatereqsc( S1, S2 )
1043	implicit none
1044	logical(lp) dgreatereqsc
1045	TYPE (real_8), INTENT (IN) :: S1
1046	real(dp), INTENT (IN) :: S2
1047
1048	dgreatereqsc=.FALSE.
1049	select case(s1%kind)
1050	case(m1,m3)
1051	IF(S1%r>=S2) dgreatereqsc=.TRUE.
1052	case(m2)
1053	IF((S1%t.SUB.'0')>=S2) dgreatereqsc=.TRUE.
1054	case default
1055	w_p=0
1056	w_p%nc=2
1057	w_p%fc='((1X,A72,/,1x,a72))'
1058	w_p%fi='(1((1X,i4)))'
1059	w_p%c(1)= " trouble in dgreatereqsc "
1060	w_p%c(2)= "s1%kind "
1061	w_p=(/s1%kind/)
1062	! call !write_e(0)
1063	end select
1064
1065	END FUNCTION dgreatereqsc
1066
1067	FUNCTION dscgreatereq( S2,S1 )
1068	implicit none
1069	logical(lp) dscgreatereq
1070	TYPE (real_8), INTENT (IN) :: S1
1071	real(dp), INTENT (IN) :: S2
1072
1073	dscgreatereq=.FALSE.
1074	select case(s1%kind)
1075	case(m1,m3)
1076	IF(S2>=S1%r) dscgreatereq=.TRUE.
1077	case(m2)
1078	IF(S2>=(S1%t.SUB.'0')) dscgreatereq=.TRUE.
1079	case default
1080
1081	w_p=0
1082	w_p%nc=2
1083	w_p%fc='((1X,A72,/,1x,a72))'
1084	w_p%fi='(1((1X,i4)))'
1085	w_p%c(1)= " trouble in dscgreatereq "
1086	w_p%c(2)= "s1%kind "
1087	w_p=(/s1%kind/)
1088	! call !write_e(0)
1089	end select
1090
1091	END FUNCTION dscgreatereq
1092
1093	FUNCTION greatereqsc( S1, S2 )
1094	implicit none
1095	logical(lp) greatereqsc
1096	TYPE (real_8), INTENT (IN) :: S1
1097	real(sp), INTENT (IN) :: S2
1098
1099	if(real_warning) call real_stop
1100	greatereqsc=.FALSE.
1101	select case(s1%kind)
1102	case(m1,m3)
1103	IF(S1%r>=S2) greatereqsc=.TRUE.
1104	case(m2)
1105	IF((S1%t.SUB.'0')>=S2) greatereqsc=.TRUE.
1106	case default
1107	w_p=0
1108	w_p%nc=2
1109	w_p%fc='((1X,A72,/,1x,a72))'
1110	w_p%fi='(1((1X,i4)))'
1111	w_p%c(1)= " trouble in greatereqsc "
1112	w_p%c(2)= "s1%kind "
1113	w_p=(/s1%kind/)
1114	! call !write_e(0)
1115	end select
1116
1117	END FUNCTION greatereqsc
1118
1119	FUNCTION scgreatereq( S2,S1 )
1120	implicit none
1121	logical(lp) scgreatereq
1122	TYPE (real_8), INTENT (IN) :: S1
1123	real(sp), INTENT (IN) :: S2
1124
1125	if(real_warning) call real_stop
1126	scgreatereq=.FALSE.
1127	select case(s1%kind)
1128	case(m1,m3)
1129	IF(S2>=S1%r) scgreatereq=.TRUE.
1130	case(m2)
1131	IF(S2>=(S1%t.SUB.'0')) scgreatereq=.TRUE.
1132	case default
1133	w_p=0
1134	w_p%nc=2
1135	w_p%fc='((1X,A72,/,1x,a72))'
1136	w_p%fi='(1((1X,i4)))'
1137	w_p%c(1)= " trouble in scgreatereq "
1138	w_p%c(2)= "s1%kind "
1139	w_p=(/s1%kind/)
1140	! call !write_e(0)
1141	end select
1142
1143	END FUNCTION scgreatereq
1144
1145	FUNCTION igreatereqsc( S1, S2 )
1146	implicit none
1147	logical(lp) igreatereqsc
1148	TYPE (real_8), INTENT (IN) :: S1
1149	integer, INTENT (IN) :: S2
1150
1151	igreatereqsc=.FALSE.
1152	select case(s1%kind)
1153	case(m1,m3)
1154	IF(S1%r>=S2) igreatereqsc=.TRUE.
1155	case(m2)
1156	IF((S1%t.SUB.'0')>=S2) igreatereqsc=.TRUE.
1157	case default
1158	w_p=0
1159	w_p%nc=2
1160	w_p%fc='((1X,A72,/,1x,a72))'
1161	w_p%fi='(1((1X,i4)))'
1162	w_p%c(1)= " trouble in igreatereqsc "
1163	w_p%c(2)= "s1%kind "
1164	w_p=(/s1%kind/)
1165	! call !write_e(0)
1166	end select
1167
1168	END FUNCTION igreatereqsc
1169
1170	FUNCTION iscgreatereq( S2,S1 )
1171	implicit none
1172	logical(lp) iscgreatereq
1173	TYPE (real_8), INTENT (IN) :: S1
1174	integer, INTENT (IN) :: S2
1175
1176	iscgreatereq=.FALSE.
1177	select case(s1%kind)
1178	case(m1,m3)
1179	IF(S2>=S1%r) iscgreatereq=.TRUE.
1180	case(m2)
1181	IF(S2>=(S1%t.SUB.'0')) iscgreatereq=.TRUE.
1182	case default
1183	w_p=0
1184	w_p%nc=2
1185	w_p%fc='((1X,A72,/,1x,a72))'
1186	w_p%fi='(1((1X,i4)))'
1187	w_p%c(1)= " trouble in iscgreatereq "
1188	w_p%c(2)= "s1%kind "
1189	w_p=(/s1%kind/)
1190	! call !write_e(0)
1191	end select
1192
1193	END FUNCTION iscgreatereq
1194
1195	FUNCTION greaterthan( S1, S2 )
1196	implicit none
1197	logical(lp) greaterthan
1198	TYPE (real_8), INTENT (IN) :: S1,S2
1199
1200	greaterthan=.FALSE.
1201	select case(s1%kind+ms*s2%kind)
1202	case(m11)
1203	IF(S1%r>S2%r) greaterthan=.TRUE.
1204	case(m12,m21,m22)
1205	select case(s1%kind+ms*s2%kind)
1206	case(m21)
1207	IF((S1%t.SUB.'0')>S2%r) greaterthan=.TRUE.
1208	case(m12)
1209	IF(S1%r>(S2%t.SUB.'0')) greaterthan=.TRUE.
1210	case(m22)
1211	IF((S1%t.SUB.'0')>(S2%t.SUB.'0')) greaterthan=.TRUE.
1212	end select
1213	case(m13,m31,m32,m23,m33)
1214	select case(s1%kind+ms*s2%kind)
1215	case(m31,m13,m33)
1216	IF(S1%r>S2%r) greaterthan=.TRUE.
1217	case(m32)
1218	IF(S1%r>(S2%t.SUB.'0')) greaterthan=.TRUE.
1219	case(m23)
1220	IF((S1%t.SUB.'0')>S2%r) greaterthan=.TRUE.
1221	end select
1222	case default
1223	w_p=0
1224	w_p%nc=2
1225	w_p%fc='((1X,A72,/,1x,a72))'
1226	w_p%fi='(2((1X,i4)))'
1227	w_p%c(1)= " trouble in greaterthan "
1228	w_p%c(2)= "s1%kind ,s2%kind "
1229	w_p=(/s1%kind ,s2%kind/)
1230	! call !write_e(0)
1231	end select
1232
1233	END FUNCTION greaterthan
1234
1235	FUNCTION igreatersc( S1, S2 )
1236	implicit none
1237	logical(lp) igreatersc
1238	TYPE (real_8), INTENT (IN) :: S1
1239	integer, INTENT (IN) :: S2
1240
1241	igreatersc=.FALSE.
1242	select case(s1%kind)
1243	case(m1,m3)
1244	IF(S1%r>S2) igreatersc=.TRUE.
1245	case(m2)
1246	IF((S1%t.SUB.'0')>S2) igreatersc=.TRUE.
1247	case default
1248	w_p=0
1249	w_p%nc=2
1250	w_p%fc='((1X,A72,/,1x,a72))'
1251	w_p%fi='(1((1X,i4)))'
1252	w_p%c(1)= " trouble in igreatersc "
1253	w_p%c(2)= "s1%kind "
1254	w_p=(/s1%kind/)
1255	! call !write_e(0)
1256	end select
1257
1258	END FUNCTION igreatersc
1259
1260	FUNCTION iscgreater( S2,S1 )
1261	implicit none
1262	logical(lp) iscgreater
1263	TYPE (real_8), INTENT (IN) :: S1
1264	integer, INTENT (IN) :: S2
1265
1266	iscgreater=.FALSE.
1267	select case(s1%kind)
1268	case(m1,m3)
1269	IF(S2>S1%r) iscgreater=.TRUE.
1270	case(m2)
1271	IF(S2>(S1%t.SUB.'0')) iscgreater=.TRUE.
1272	case default
1273	w_p=0
1274	w_p%nc=2
1275	w_p%fc='((1X,A72,/,1x,a72))'
1276	w_p%fi='(1((1X,i4)))'
1277	w_p%c(1)= " trouble in iscgreater "
1278	w_p%c(2)= "s1%kind "
1279	w_p=(/s1%kind/)
1280	! call !write_e(0)
1281	end select
1282
1283	END FUNCTION iscgreater
1284
1285	FUNCTION dgreatersc( S1, S2 )
1286	implicit none
1287	logical(lp) dgreatersc
1288	TYPE (real_8), INTENT (IN) :: S1
1289	real(dp), INTENT (IN) :: S2
1290
1291	dgreatersc=.FALSE.
1292	select case(s1%kind)
1293	case(m1,m3)
1294	IF(S1%r>S2) dgreatersc=.TRUE.
1295	case(m2)
1296	IF((S1%t.SUB.'0')>S2) dgreatersc=.TRUE.
1297	case default
1298	w_p=0
1299	w_p%nc=2
1300	w_p%fc='((1X,A72,/,1x,a72))'
1301	w_p%fi='(1((1X,i4)))'
1302	w_p%c(1)= " trouble in dgreatersc "
1303	w_p%c(2)= "s1%kind "
1304	w_p=(/s1%kind/)
1305	! call !write_e(0)
1306	end select
1307
1308	END FUNCTION dgreatersc
1309
1310	FUNCTION dscgreater( S2,S1 )
1311	implicit none
1312	integer ipause, mypause
1313	logical(lp) dscgreater
1314	TYPE (real_8), INTENT (IN) :: S1
1315	real(dp), INTENT (IN) :: S2
1316
1317	dscgreater=.FALSE.
1318	select case(s1%kind)
1319	case(m1,m3)
1320	IF(S2>S1%r) dscgreater=.TRUE.
1321	case(m2)
1322	IF(S2>(S1%t.SUB.'0')) dscgreater=.TRUE.
1323	case default
1324	w_p=0
1325	w_p%nc=2
1326	w_p%fc='((1X,A72,/,1x,a72))'
1327	w_p%fi='(1((1X,i4)))'
1328	w_p%c(1)= " trouble in dscgreater "
1329	w_p%c(2)= "s1%kind "
1330	w_p=(/s1%kind/)
1331	! call !write_e(0)
1332	ipause=mypause(0)
1333	end select
1334
1335	END FUNCTION dscgreater
1336
1337	FUNCTION greatersc( S1, S2 )
1338	implicit none
1339	logical(lp) greatersc
1340	TYPE (real_8), INTENT (IN) :: S1
1341	real(sp), INTENT (IN) :: S2
1342
1343	if(real_warning) call real_stop
1344	greatersc=.FALSE.
1345	select case(s1%kind)
1346	case(m1,m3)
1347	IF(S1%r>S2) greatersc=.TRUE.
1348	case(m2)
1349	IF((S1%t.SUB.'0')>S2) greatersc=.TRUE.
1350	case default
1351	w_p=0
1352	w_p%nc=2
1353	w_p%fc='((1X,A72,/,1x,a72))'
1354	w_p%fi='(2((1X,i4)))'
1355	w_p%c(1)= " trouble in greatersc "
1356	w_p%c(2)= "s1%kind "
1357	w_p=(/s1%kind /)
1358	! call !write_e(0)
1359	end select
1360
1361	END FUNCTION greatersc
1362
1363	FUNCTION scgreater( S2,S1 )
1364	implicit none
1365	logical(lp) scgreater
1366	TYPE (real_8), INTENT (IN) :: S1
1367	real(sp), INTENT (IN) :: S2
1368
1369	if(real_warning) call real_stop
1370	scgreater=.FALSE.
1371	select case(s1%kind)
1372	case(m1,m3)
1373	IF(S2>S1%r) scgreater=.TRUE.
1374	case(m2)
1375	IF(S2>(S1%t.SUB.'0')) scgreater=.TRUE.
1376	case default
1377	w_p=0
1378	w_p%nc=2
1379	w_p%fc='((1X,A72,/,1x,a72))'
1380	w_p%fi='(2((1X,i4)))'
1381	w_p%c(1)= " trouble in scgreater "
1382	w_p%c(2)= "s1%kind "
1383	w_p=(/s1%kind /)
1384	! call !write_e(0)
1385	end select
1386
1387	END FUNCTION scgreater
1388
1389	FUNCTION lessthan( S1, S2 )
1390	implicit none
1391	logical(lp) lessthan
1392	TYPE (real_8), INTENT (IN) :: S1,S2
1393
1394	lessthan=.FALSE.
1395	select case(s1%kind+ms*s2%kind)
1396	case(m11)
1397	IF(S1%r<S2%r) lessthan=.TRUE.
1398	case(m12,m21,m22)
1399	select case(s1%kind+ms*s2%kind)
1400	case(m21)
1401	IF((S1%t.SUB.'0')<S2%r) lessthan=.TRUE.
1402	case(m12)
1403	IF(S1%r<(S2%t.SUB.'0')) lessthan=.TRUE.
1404	case(m22)
1405	IF((S1%t.SUB.'0')<(S2%t.SUB.'0')) lessthan=.TRUE.
1406	end select
1407	case(m13,m31,m32,m23,m33)
1408	select case(s1%kind+ms*s2%kind)
1409	case(m31,m13,m33)
1410	IF(S1%r<S2%r) lessthan=.TRUE.
1411	case(m32)
1412	IF(S1%r<(S2%t.SUB.'0')) lessthan=.TRUE.
1413	case(m23)
1414	IF((S1%t.SUB.'0')<S2%r) lessthan=.TRUE.
1415	end select
1416	case default
1417	w_p=0
1418	w_p%nc=2
1419	w_p%fc='((1X,A72,/,1x,a72))'
1420	w_p%fi='(2((1X,i4)))'
1421	w_p%c(1)= " trouble in lessthan "
1422	w_p%c(2)= "s1%kind ,s2%kind "
1423	w_p=(/s1%kind ,s2%kind/)
1424	! call !write_e(0)
1425	end select
1426
1427	END FUNCTION lessthan
1428
1429	FUNCTION dlessthansc( S1, S2 )
1430	implicit none
1431	logical(lp) dlessthansc
1432	TYPE (real_8), INTENT (IN) :: S1
1433	real(dp), INTENT (IN) :: S2
1434
1435	dlessthansc=.FALSE.
1436	select case(s1%kind)
1437	case(m1,m3)
1438	IF(S1%r<S2) dlessthansc=.TRUE.
1439	case(m2)
1440	IF((S1%t.SUB.'0')<S2) dlessthansc=.TRUE.
1441	case default
1442	w_p=0
1443	w_p%nc=2
1444	w_p%fc='((1X,A72,/,1x,a72))'
1445	w_p%fi='(2((1X,i4)))'
1446	w_p%c(1)= " trouble in dlessthansc "
1447	w_p%c(2)= "s1%kind "
1448	w_p=(/s1%kind /)
1449	! call !write_e(0)
1450	end select
1451
1452	END FUNCTION dlessthansc
1453
1454
1455	FUNCTION dsclessthan( S2,S1 )
1456	implicit none
1457	logical(lp) dsclessthan
1458	TYPE (real_8), INTENT (IN) :: S1
1459	real(dp), INTENT (IN) :: S2
1460
1461	dsclessthan=.FALSE.
1462	select case(s1%kind)
1463	case(m1,m3)
1464	IF(S2<S1%r) dsclessthan=.TRUE.
1465	case(m2)
1466	IF(S2<(S1%t.SUB.'0')) dsclessthan=.TRUE.
1467	case default
1468	w_p=0
1469	w_p%nc=2
1470	w_p%fc='((1X,A72,/,1x,a72))'
1471	w_p%fi='(2((1X,i4)))'
1472	w_p%c(1)= " trouble in dsclessthan "
1473	w_p%c(2)= "s1%kind "
1474	w_p=(/s1%kind /)
1475	! call !write_e(0)
1476	end select
1477
1478	END FUNCTION dsclessthan
1479
1480	FUNCTION lessthansc( S1, S2 )
1481	implicit none
1482	logical(lp) lessthansc
1483	TYPE (real_8), INTENT (IN) :: S1
1484	real(sp), INTENT (IN) :: S2
1485
1486	if(real_warning) call real_stop
1487	lessthansc=.FALSE.
1488	select case(s1%kind)
1489	case(m1,m3)
1490	IF(S1%r<S2) lessthansc=.TRUE.
1491	case(m2)
1492	IF((S1%t.SUB.'0')<S2) lessthansc=.TRUE.
1493	case default
1494	w_p=0
1495	w_p%nc=2
1496	w_p%fc='((1X,A72,/,1x,a72))'
1497	w_p%fi='(2((1X,i4)))'
1498	w_p%c(1)= " trouble in lessthansc "
1499	w_p%c(2)= "s1%kind "
1500	w_p=(/s1%kind /)
1501	! call !write_e(0)
1502	end select
1503
1504	END FUNCTION lessthansc
1505
1506
1507	FUNCTION sclessthan( S2,S1 )
1508	implicit none
1509	logical(lp) sclessthan
1510	TYPE (real_8), INTENT (IN) :: S1
1511	real(sp), INTENT (IN) :: S2
1512
1513	if(real_warning) call real_stop
1514	sclessthan=.FALSE.
1515	select case(s1%kind)
1516	case(m1,m3)
1517	IF(S2<S1%r) sclessthan=.TRUE.
1518	case(m2)
1519	IF(S2<(S1%t.SUB.'0')) sclessthan=.TRUE.
1520	case default
1521	w_p=0
1522	w_p%nc=2
1523	w_p%fc='((1X,A72,/,1x,a72))'
1524	w_p%fi='(2((1X,i4)))'
1525	w_p%c(1)= " trouble in sclessthan "
1526	w_p%c(2)= "s1%kind "
1527	w_p=(/s1%kind /)
1528	! call !write_e(0)
1529	end select
1530
1531	END FUNCTION sclessthan
1532
1533	FUNCTION ilessthansc( S1, S2 )
1534	implicit none
1535	logical(lp) ilessthansc
1536	TYPE (real_8), INTENT (IN) :: S1
1537	integer, INTENT (IN) :: S2
1538
1539	ilessthansc=.FALSE.
1540	select case(s1%kind)
1541	case(m1,m3)
1542	IF(S1%r<S2) ilessthansc=.TRUE.
1543	case(m2)
1544	IF((S1%t.SUB.'0')<S2) ilessthansc=.TRUE.
1545	case default
1546	w_p=0
1547	w_p%nc=2
1548	w_p%fc='((1X,A72,/,1x,a72))'
1549	w_p%fi='(2((1X,i4)))'
1550	w_p%c(1)= " trouble in ilessthansc "
1551	w_p%c(2)= "s1%kind ,s2%kind "
1552	w_p%c(2)= "s1%kind "
1553	w_p=(/s1%kind /)
1554	end select
1555
1556	END FUNCTION ilessthansc
1557
1558
1559	FUNCTION isclessthan( S2,S1 )
1560	implicit none
1561	logical(lp) isclessthan
1562	TYPE (real_8), INTENT (IN) :: S1
1563	integer, INTENT (IN) :: S2
1564
1565	isclessthan=.FALSE.
1566	select case(s1%kind)
1567	case(m1,m3)
1568	IF(S2<S1%r) isclessthan=.TRUE.
1569	case(m2)
1570	IF(S2<(S1%t.SUB.'0')) isclessthan=.TRUE.
1571	case default
1572	w_p=0
1573	w_p%nc=2
1574	w_p%fc='((1X,A72,/,1x,a72))'
1575	w_p%fi='(2((1X,i4)))'
1576	w_p%c(1)= " trouble in isclessthan "
1577	w_p%c(2)= "s1%kind "
1578	w_p=(/s1%kind /)
1579	! call !write_e(0)
1580	end select
1581
1582	END FUNCTION isclessthan
1583
1584	FUNCTION lesseq( S1, S2 )
1585	implicit none
1586	logical(lp) lesseq
1587	TYPE (real_8), INTENT (IN) :: S1,S2
1588
1589	lesseq=.FALSE.
1590	select case(s1%kind+ms*s2%kind)
1591	case(m11)
1592	IF(S1%r<=S2%r) lesseq=.TRUE.
1593	case(m12,m21,m22)
1594	select case(s1%kind+ms*s2%kind)
1595	case(m21)
1596	IF((S1%t.SUB.'0')<=S2%r) lesseq=.TRUE.
1597	case(m12)
1598	IF(S1%r<=(S2%t.SUB.'0')) lesseq=.TRUE.
1599	case(m22)
1600	IF((S1%t.SUB.'0')<=(S2%t.SUB.'0')) lesseq=.TRUE.
1601	end select
1602	case(m13,m31,m32,m23,m33)
1603	select case(s1%kind+ms*s2%kind)
1604	case(m31,m13,m33)
1605	IF(S1%r<=S2%r) lesseq=.TRUE.
1606	case(m32)
1607	IF(S1%r<=(S2%t.SUB.'0')) lesseq=.TRUE.
1608	case(m23)
1609	IF((S1%t.SUB.'0')<=S2%r) lesseq=.TRUE.
1610	end select
1611	case default
1612	w_p=0
1613	w_p%nc=2
1614	w_p%fc='((1X,A72,/,1x,a72))'
1615	w_p%fi='(2((1X,i4)))'
1616	w_p%c(1)= " trouble in lesseq "
1617	w_p%c(2)= "s1%kind ,s2%kind "
1618	w_p=(/s1%kind ,s2%kind/)
1619	! call !write_e(0)
1620	end select
1621
1622
1623	END FUNCTION lesseq
1624
1625	FUNCTION dlesseqsc( S1, S2 )
1626	implicit none
1627	logical(lp) dlesseqsc
1628	TYPE (real_8), INTENT (IN) :: S1
1629	real(dp), INTENT (IN) :: S2
1630
1631	dlesseqsc=.FALSE.
1632	select case(s1%kind)
1633	case(m1,m3)
1634	IF(S1%r<=S2) dlesseqsc=.TRUE.
1635	case(m2)
1636	IF((S1%t.SUB.'0')<=S2) dlesseqsc=.TRUE.
1637	case default
1638	w_p=0
1639	w_p%nc=2
1640	w_p%fc='((1X,A72,/,1x,a72))'
1641	w_p%fi='(2((1X,i4)))'
1642	w_p%c(1)= " trouble in dlesseqsc "
1643	w_p%c(2)= "s1%kind "
1644	w_p=(/s1%kind /)
1645	! call !write_e(0)
1646	end select
1647
1648	END FUNCTION dlesseqsc
1649
1650	FUNCTION dsclesseq( S2,S1 )
1651	implicit none
1652	logical(lp) dsclesseq
1653	TYPE (real_8), INTENT (IN) :: S1
1654	real(dp), INTENT (IN) :: S2
1655
1656	dsclesseq=.FALSE.
1657	select case(s1%kind)
1658	case(m1,m3)
1659	IF(S2<=S1%r) dsclesseq=.TRUE.
1660	case(m2)
1661	IF(S2<=(S1%t.SUB.'0')) dsclesseq=.TRUE.
1662	case default
1663	w_p=0
1664	w_p%nc=2
1665	w_p%fc='((1X,A72,/,1x,a72))'
1666	w_p%fi='(2((1X,i4)))'
1667	w_p%c(1)= " trouble in dsclesseq "
1668	w_p%c(2)= "s1%kind "
1669	w_p=(/s1%kind /)
1670	! call !write_e(0)
1671	end select
1672
1673	END FUNCTION dsclesseq
1674
1675	FUNCTION lesseqsc( S1, S2 )
1676	implicit none
1677	logical(lp) lesseqsc
1678	TYPE (real_8), INTENT (IN) :: S1
1679	real(sp), INTENT (IN) :: S2
1680
1681	if(real_warning) call real_stop
1682	lesseqsc=.FALSE.
1683	select case(s1%kind)
1684	case(m1,m3)
1685	IF(S1%r<=S2) lesseqsc=.TRUE.
1686	case(m2)
1687	IF((S1%t.SUB.'0')<=S2) lesseqsc=.TRUE.
1688	case default
1689	w_p=0
1690	w_p%nc=2
1691	w_p%fc='((1X,A72,/,1x,a72))'
1692	w_p%fi='(2((1X,i4)))'
1693	w_p%c(1)= " trouble in lesseqsc "
1694	w_p%c(2)= "s1%kind "
1695	w_p=(/s1%kind /)
1696	! call !write_e(0)
1697	end select
1698
1699	END FUNCTION lesseqsc
1700
1701	FUNCTION sclesseq( S2,S1 )
1702	implicit none
1703	logical(lp) sclesseq
1704	TYPE (real_8), INTENT (IN) :: S1
1705	real(sp), INTENT (IN) :: S2
1706
1707	if(real_warning) call real_stop
1708	sclesseq=.FALSE.
1709	select case(s1%kind)
1710	case(m1,m3)
1711	IF(S2<=S1%r) sclesseq=.TRUE.
1712	case(m2)
1713	IF(S2<=(S1%t.SUB.'0')) sclesseq=.TRUE.
1714	case default
1715	w_p=0
1716	w_p%nc=2
1717	w_p%fc='((1X,A72,/,1x,a72))'
1718	w_p%fi='(2((1X,i4)))'
1719	w_p%c(1)= " trouble in sclesseq "
1720	w_p%c(2)= "s1%kind "
1721	w_p=(/s1%kind /)
1722	! call !write_e(0)
1723	end select
1724
1725	END FUNCTION sclesseq
1726
1727	FUNCTION ilesseqsc( S1, S2 )
1728	implicit none
1729	logical(lp) ilesseqsc
1730	TYPE (real_8), INTENT (IN) :: S1
1731	integer, INTENT (IN) :: S2
1732
1733	ilesseqsc=.FALSE.
1734	select case(s1%kind)
1735	case(m1,m3)
1736	IF(S1%r<=S2) ilesseqsc=.TRUE.
1737	case(m2)
1738	IF((S1%t.SUB.'0')<=S2) ilesseqsc=.TRUE.
1739	case default
1740	w_p=0
1741	w_p%nc=2
1742	w_p%fc='((1X,A72,/,1x,a72))'
1743	w_p%fi='(2((1X,i4)))'
1744	w_p%c(1)= " trouble in ilesseqsc "
1745	w_p%c(2)= "s1%kind "
1746	w_p=(/s1%kind /)
1747	! call !write_e(0)
1748	end select
1749
1750	END FUNCTION ilesseqsc
1751
1752	FUNCTION isclesseq( S2,S1 )
1753	implicit none
1754	logical(lp) isclesseq
1755	TYPE (real_8), INTENT (IN) :: S1
1756	integer, INTENT (IN) :: S2
1757
1758	isclesseq=.FALSE.
1759	select case(s1%kind)
1760	case(m1,m3)
1761	IF(S2<=S1%r) isclesseq=.TRUE.
1762	case(m2)
1763	IF(S2<=(S1%t.SUB.'0')) isclesseq=.TRUE.
1764	case default
1765	w_p=0
1766	w_p%nc=2
1767	w_p%fc='((1X,A72,/,1x,a72))'
1768	w_p%fi='(2((1X,i4)))'
1769	w_p%c(1)= " trouble in isclesseq "
1770	w_p%c(2)= "s1%kind "
1771	w_p=(/s1%kind /)
1772	! call !write_e(0)
1773	end select
1774
1775	END FUNCTION isclesseq
1776
1777	FUNCTION eq( S1, S2 )
1778	implicit none
1779	logical(lp) eq
1780	TYPE (real_8), INTENT (IN) :: S1,S2
1781
1782	eq=.FALSE.
1783	select case(s1%kind+ms*s2%kind)
1784	case(m11)
1785	IF(S1%r==S2%r) eq=.TRUE.
1786	case(m12,m21,m22)
1787	select case(s1%kind+ms*s2%kind)
1788	case(m21)
1789	IF((S1%t.SUB.'0')==S2%r) eq=.TRUE.
1790	case(m12)
1791	IF(S1%r==(S2%t.SUB.'0')) eq=.TRUE.
1792	case(m22)
1793	IF((S1%t.SUB.'0')==(S2%t.SUB.'0')) eq=.TRUE.
1794	end select
1795	case(m13,m31,m32,m23,m33)
1796	select case(s1%kind+ms*s2%kind)
1797	case(m31,m13,m33)
1798	IF(S1%r==S2%r) eq=.TRUE.
1799	case(m32)
1800	IF(S1%r==(S2%t.SUB.'0')) eq=.TRUE.
1801	case(m23)
1802	IF((S1%t.SUB.'0')==S2%r) eq=.TRUE.
1803	end select
1804	case default
1805	w_p=0
1806	w_p%nc=2
1807	w_p%fc='((1X,A72,/,1x,a72))'
1808	w_p%fi='(2((1X,i4)))'
1809	w_p%c(1)= " trouble in eq "
1810	w_p%c(2)= "s1%kind ,s2%kind "
1811	w_p=(/s1%kind ,s2%kind/)
1812	! call !write_e(0)
1813	end select
1814
1815
1816	END FUNCTION eq
1817
1818	FUNCTION deqsc( S1, S2 )
1819	implicit none
1820	logical(lp) deqsc
1821	TYPE (real_8), INTENT (IN) :: S1
1822	real(dp), INTENT (IN) :: S2
1823
1824	deqsc=.FALSE.
1825	select case(s1%kind)
1826	case(m1,m3)
1827	IF(S1%r==S2) deqsc=.TRUE.
1828	case(m2)
1829	IF((S1%t.SUB.'0')==S2) deqsc=.TRUE.
1830	case default
1831	w_p=0
1832	w_p%nc=2
1833	w_p%fc='((1X,A72,/,1x,a72))'
1834	w_p%fi='(2((1X,i4)))'
1835	w_p%c(1)= " trouble in deqsc "
1836	w_p%c(2)= "s1%kind "
1837	w_p=(/s1%kind /)
1838	! call !write_e(0)
1839	end select
1840
1841	END FUNCTION deqsc
1842
1843	FUNCTION dsceq( S2,S1 )
1844	implicit none
1845	logical(lp) dsceq
1846	TYPE (real_8), INTENT (IN) :: S1
1847	real(dp), INTENT (IN) :: S2
1848
1849	dsceq=.FALSE.
1850	select case(s1%kind)
1851	case(m1,m3)
1852	IF(S2==S1%r) dsceq=.TRUE.
1853	case(m2)
1854	IF(S2==(S1%t.SUB.'0')) dsceq=.TRUE.
1855	case default
1856	w_p=0
1857	w_p%nc=2
1858	w_p%fc='((1X,A72,/,1x,a72))'
1859	w_p%fi='(2((1X,i4)))'
1860	w_p%c(1)= " trouble in dsceq "
1861	w_p%c(2)= "s1%kind "
1862	w_p=(/s1%kind /)
1863	! call !write_e(0)
1864	end select
1865
1866	END FUNCTION dsceq
1867
1868	FUNCTION eqsc( S1, S2 )
1869	implicit none
1870	logical(lp) eqsc
1871	TYPE (real_8), INTENT (IN) :: S1
1872	real(sp), INTENT (IN) :: S2
1873
1874	if(real_warning) call real_stop
1875	eqsc=.FALSE.
1876	select case(s1%kind)
1877	case(m1,m3)
1878	IF(S1%r==S2) eqsc=.TRUE.
1879	case(m2)
1880	IF((S1%t.SUB.'0')==S2) eqsc=.TRUE.
1881	case default
1882	w_p=0
1883	w_p%nc=2
1884	w_p%fc='((1X,A72,/,1x,a72))'
1885	w_p%fi='(2((1X,i4)))'
1886	w_p%c(1)= " trouble in eqsc "
1887	w_p%c(2)= "s1%kind "
1888	w_p=(/s1%kind /)
1889	! call !write_e(0)
1890	end select
1891
1892	END FUNCTION eqsc
1893
1894	FUNCTION sceq( S2,S1 )
1895	implicit none
1896	logical(lp) sceq
1897	TYPE (real_8), INTENT (IN) :: S1
1898	real(sp), INTENT (IN) :: S2
1899
1900	if(real_warning) call real_stop
1901	sceq=.FALSE.
1902	select case(s1%kind)
1903	case(m1,m3)
1904	IF(S2==S1%r) sceq=.TRUE.
1905	case(m2)
1906	IF(S2==(S1%t.SUB.'0')) sceq=.TRUE.
1907	case default
1908	w_p=0
1909	w_p%nc=2
1910	w_p%fc='((1X,A72,/,1x,a72))'
1911	w_p%fi='(2((1X,i4)))'
1912	w_p%c(1)= " trouble in sceq "
1913	w_p%c(2)= "s1%kind "
1914	w_p=(/s1%kind /)
1915	! call !write_e(0)
1916	end select
1917
1918	END FUNCTION sceq
1919
1920
1921	FUNCTION ieqsc( S1, S2 )
1922	implicit none
1923	logical(lp) ieqsc
1924	TYPE (real_8), INTENT (IN) :: S1
1925	integer, INTENT (IN) :: S2
1926
1927	ieqsc=.FALSE.
1928	select case(s1%kind)
1929	case(m1,m3)
1930	IF(S1%r==S2) ieqsc=.TRUE.
1931	case(m2)
1932	IF((S1%t.SUB.'0')==S2) ieqsc=.TRUE.
1933	case default
1934	w_p=0
1935	w_p%nc=2
1936	w_p%fc='((1X,A72,/,1x,a72))'
1937	w_p%fi='(2((1X,i4)))'
1938	w_p%c(1)= " trouble in ieqsc "
1939	w_p%c(2)= "s1%kind "
1940	w_p=(/s1%kind /)
1941	! call !write_e(0)
1942	end select
1943
1944	END FUNCTION ieqsc
1945
1946	FUNCTION isceq( S2,S1 )
1947	implicit none
1948	logical(lp) isceq
1949	TYPE (real_8), INTENT (IN) :: S1
1950	integer, INTENT (IN) :: S2
1951
1952	isceq=.FALSE.
1953	select case(s1%kind)
1954	case(m1,m3)
1955	IF(S2==S1%r) isceq=.TRUE.
1956	case(m2)
1957	IF(S2==(S1%t.SUB.'0')) isceq=.TRUE.
1958	case default
1959	w_p=0
1960	w_p%nc=2
1961	w_p%fc='((1X,A72,/,1x,a72))'
1962	w_p%fi='(2((1X,i4)))'
1963	w_p%c(1)= " trouble in isceq "
1964	w_p%c(2)= "s1%kind "
1965	w_p=(/s1%kind /)
1966	! call !write_e(0)
1967	end select
1968
1969	END FUNCTION isceq
1970
1971	FUNCTION neq( S1, S2 )
1972	implicit none
1973	logical(lp) neq
1974	TYPE (real_8), INTENT (IN) :: S1,S2
1975
1976	neq=.FALSE.
1977	select case(s1%kind+ms*s2%kind)
1978	case(m11)
1979	IF(S1%r/=S2%r) neq=.TRUE.
1980	case(m12,m21,m22)
1981	select case(s1%kind+ms*s2%kind)
1982	case(m21)
1983	IF((S1%t.SUB.'0')/=S2%r) neq=.TRUE.
1984	case(m12)
1985	IF(S1%r/=(S2%t.SUB.'0')) neq=.TRUE.
1986	case(m22)
1987	IF((S1%t.SUB.'0')/=(S2%t.SUB.'0')) neq=.TRUE.
1988	end select
1989	case(m13,m31,m32,m23,m33)
1990	select case(s1%kind+ms*s2%kind)
1991	case(m31,m13,m33)
1992	IF(S1%r/=S2%r) neq=.TRUE.
1993	case(m32)
1994	IF(S1%r/=(S2%t.SUB.'0')) neq=.TRUE.
1995	case(m23)
1996	IF((S1%t.SUB.'0')/=S2%r) neq=.TRUE.
1997	end select
1998	case default
1999	w_p=0
2000	w_p%nc=2
2001	w_p%fc='((1X,A72,/,1x,a72))'
2002	w_p%fi='(2((1X,i4)))'
2003	w_p%c(1)= " trouble in neq "
2004	w_p%c(2)= "s1%kind ,s2%kind "
2005	w_p=(/s1%kind ,s2%kind/)
2006	! call !write_e(0)
2007	end select
2008
2009
2010	END FUNCTION neq
2011
2012	FUNCTION dneqsc( S1, S2 )
2013	implicit none
2014	logical(lp) dneqsc
2015	TYPE (real_8), INTENT (IN) :: S1
2016	real(dp), INTENT (IN) :: S2
2017
2018	dneqsc=.FALSE.
2019	select case(s1%kind)
2020	case(m1,m3)
2021	IF(S1%r/=S2) dneqsc=.TRUE.
2022	case(m2)
2023	IF((S1%t.SUB.'0')/=S2) dneqsc=.TRUE.
2024	case default
2025	w_p=0
2026	w_p%nc=2
2027	w_p%fc='((1X,A72,/,1x,a72))'
2028	w_p%fi='(2((1X,i4)))'
2029	w_p%c(1)= " trouble in dneqsc "
2030	w_p%c(2)= "s1%kind "
2031	w_p=(/s1%kind /)
2032	! call !write_e(0)
2033	end select
2034
2035	END FUNCTION dneqsc
2036
2037	FUNCTION dscneq( S2,S1 )
2038	implicit none
2039	logical(lp) dscneq
2040	TYPE (real_8), INTENT (IN) :: S1
2041	real(dp), INTENT (IN) :: S2
2042
2043	dscneq=.FALSE.
2044	select case(s1%kind)
2045	case(m1,m3)
2046	IF(S2/=S1%r) dscneq=.TRUE.
2047	case(m2)
2048	IF(S2/=(S1%t.SUB.'0')) dscneq=.TRUE.
2049	case default
2050	w_p=0
2051	w_p%nc=2
2052	w_p%fc='((1X,A72,/,1x,a72))'
2053	w_p%fi='(2((1X,i4)))'
2054	w_p%c(1)= " trouble in dscneq "
2055	w_p%c(2)= "s1%kind "
2056	w_p=(/s1%kind /)
2057	! call !write_e(0)
2058	end select
2059
2060	END FUNCTION dscneq
2061
2062	FUNCTION neqsc( S1, S2 )
2063	implicit none
2064	logical(lp) neqsc
2065	TYPE (real_8), INTENT (IN) :: S1
2066	real(sp), INTENT (IN) :: S2
2067
2068	if(real_warning) call real_stop
2069	neqsc=.FALSE.
2070	select case(s1%kind)
2071	case(m1,m3)
2072	IF(S1%r/=S2) neqsc=.TRUE.
2073	case(m2)
2074	IF((S1%t.SUB.'0')/=S2) neqsc=.TRUE.
2075	case default
2076	w_p=0
2077	w_p%nc=2
2078	w_p%fc='((1X,A72,/,1x,a72))'
2079	w_p%fi='(2((1X,i4)))'
2080	w_p%c(1)= " trouble in neqsc "
2081	w_p%c(2)= "s1%kind "
2082	w_p=(/s1%kind /)
2083	! call !write_e(0)
2084	end select
2085
2086	END FUNCTION neqsc
2087
2088	FUNCTION scneq( S2,S1 )
2089	implicit none
2090	logical(lp) scneq
2091	TYPE (real_8), INTENT (IN) :: S1
2092	real(sp), INTENT (IN) :: S2
2093
2094	if(real_warning) call real_stop
2095	scneq=.FALSE.
2096	select case(s1%kind)
2097	case(m1,m3)
2098	IF(S2/=S1%r) scneq=.TRUE.
2099	case(m2)
2100	IF(S2/=(S1%t.SUB.'0')) scneq=.TRUE.
2101	case default
2102	w_p=0
2103	w_p%nc=2
2104	w_p%fc='((1X,A72,/,1x,a72))'
2105	w_p%fi='(2((1X,i4)))'
2106	w_p%c(1)= " trouble in scneq "
2107	w_p%c(2)= "s1%kind "
2108	w_p=(/s1%kind /)
2109	! call !write_e(0)
2110	end select
2111
2112	END FUNCTION scneq
2113
2114	FUNCTION ineqsc( S1, S2 )
2115	implicit none
2116	logical(lp) ineqsc
2117	TYPE (real_8), INTENT (IN) :: S1
2118	integer, INTENT (IN) :: S2
2119
2120	ineqsc=.FALSE.
2121	select case(s1%kind)
2122	case(m1,m3)
2123	IF(S1%r/=S2) ineqsc=.TRUE.
2124	case(m2)
2125	IF((S1%t.SUB.'0')/=S2) ineqsc=.TRUE.
2126	case default
2127	w_p=0
2128	w_p%nc=2
2129	w_p%fc='((1X,A72,/,1x,a72))'
2130	w_p%fi='(2((1X,i4)))'
2131	w_p%c(1)= " trouble in ineqsc "
2132	w_p%c(2)= "s1%kind "
2133	w_p=(/s1%kind /)
2134	! call !write_e(0)
2135	end select
2136
2137	END FUNCTION ineqsc
2138
2139	FUNCTION iscneq( S2,S1 )
2140	implicit none
2141	logical(lp) iscneq
2142	TYPE (real_8), INTENT (IN) :: S1
2143	integer, INTENT (IN) :: S2
2144
2145	iscneq=.FALSE.
2146	select case(s1%kind)
2147	case(m1,m3)
2148	IF(S2/=S1%r) iscneq=.TRUE.
2149	case(m2)
2150	IF(S2/=(S1%t.SUB.'0')) iscneq=.TRUE.
2151	case default
2152	w_p=0
2153	w_p%nc=2
2154	w_p%fc='((1X,A72,/,1x,a72))'
2155	w_p%fi='(2((1X,i4)))'
2156	w_p%c(1)= " trouble in iscneq "
2157	w_p%c(2)= "s1%kind "
2158	w_p=(/s1%kind /)
2159	! call !write_e(0)
2160	end select
2161
2162	END FUNCTION iscneq
2163
2164	FUNCTION dexpt( S1 )
2165	implicit none
2166	TYPE (real_8) dexpt
2167	TYPE (real_8), INTENT (IN) :: S1
2168	integer localmaster
2169
2170	select case(s1%kind)
2171	case(m1)
2172	dexpt%r=exp(s1%r)
2173	dexpt%kind=1
2174	case(m2)
2175	localmaster=master
2176	call ass(dexpt)
2177	dexpt%t= exp(s1%t)
2178	master=localmaster
2179	case(m3)
2180	if(knob) then
2181	localmaster=master
2182	call ass(dexpt)
2183
2184	call varfk1(S1)
2185	dexpt%t= exp(varf1)
2186	master=localmaster
2187	else
2188	dexpt%r= exp(S1%r)
2189	dexpt%kind=1
2190	endif
2191
2192	case default
2193	w_p=0
2194	w_p%nc=2
2195	w_p%fc='((1X,A72,/,1x,a72))'
2196	w_p%fi='(2((1X,i4)))'
2197	w_p%c(1)= " trouble in dexpt "
2198	w_p%c(2)= "s1%kind "
2199	w_p=(/s1%kind /)
2200	! call !write_e(0)
2201	end select
2202	END FUNCTION dexpt
2203
2204	FUNCTION abst( S1 )
2205	implicit none
2206	real(dp) abst
2207	TYPE (real_8), INTENT (IN) :: S1
2208
2209
2210	select case(s1%kind)
2211	case(m1,m3)
2212	abst=abs(s1%r)
2213	case(m2)
2214	abst=abs(s1%t.sub.'0')
2215	case default
2216	w_p=0
2217	w_p%nc=2
2218	w_p%fc='((1X,A72,/,1x,a72))'
2219	w_p%fi='(2((1X,i4)))'
2220	w_p%c(1)= " trouble in abst "
2221	w_p%c(2)= "s1%kind "
2222	w_p=(/s1%kind /)
2223	! call !write_e(0)
2224	end select
2225	END FUNCTION abst
2226
2227	FUNCTION Pabs( S1 )
2228	implicit none
2229	TYPE (real_8) Pabs
2230	TYPE (real_8), INTENT (IN) :: S1
2231	integer localmaster
2232
2233	select case(s1%kind)
2234	case(m1)
2235	Pabs%R=abs(s1%r)
2236	Pabs%kind=1
2237	case(m2)
2238	localmaster=master
2239	call ass(Pabs)
2240	IF((s1%t.sub.'0')<0) THEN
2241	Pabs%T=-s1%t
2242	ELSE
2243	Pabs%T=s1%t
2244	ENDIF
2245	master=localmaster
2246	case(m3)
2247	if(knob) then
2248	localmaster=master
2249	call ass(Pabs)
2250	call varfk1(S1)
2251	IF((s1%t.sub.'0')<0) THEN
2252	Pabs%T=-varf1
2253	ELSE
2254	Pabs%T=varf1
2255	ENDIF
2256	master=localmaster
2257	else
2258	PABS%r= ABS(S1%r)
2259	PABS%kind=1
2260	endif
2261	case default
2262	w_p=0
2263	w_p%nc=2
2264	w_p%fc='((1X,A72,/,1x,a72))'
2265	w_p%fi='(2((1X,i4)))'
2266	w_p%c(1)= " trouble in Pabs "
2267	w_p%c(2)= "s1%kind "
2268	w_p=(/s1%kind /)
2269	! call !write_e(0)
2270	end select
2271	END FUNCTION Pabs
2272
2273
2274
2275	FUNCTION full_abst( S1 )
2276	implicit none
2277	real(dp) full_abst
2278	TYPE (real_8), INTENT (IN) :: S1
2279
2280
2281	select case(s1%kind)
2282	case(m1,m3)
2283	full_abst=abs(s1%r)
2284	case(m2)
2285	full_abst=full_abs(s1%t) ! 2002.10.17
2286	case default
2287	w_p=0
2288	w_p%nc=2
2289	w_p%fc='((1X,A72,/,1x,a72))'
2290	w_p%fi='(2((1X,i4)))'
2291	w_p%c(1)= " trouble in full_abst "
2292	w_p%c(2)= "s1%kind "
2293	w_p=(/s1%kind /)
2294	! call !write_e(0)
2295	end select
2296	END FUNCTION full_abst
2297
2298	FUNCTION dtant( S1 )
2299	implicit none
2300	TYPE (real_8) dtant
2301	TYPE (real_8), INTENT (IN) :: S1
2302	integer localmaster
2303
2304	select case(s1%kind)
2305	case(m1)
2306	dtant%r=tan(s1%r)
2307	dtant%kind=1
2308	case(m2)
2309	localmaster=master
2310	call ass(dtant)
2311	dtant%t= tan(s1%t)
2312	master=localmaster
2313	case(m3)
2314	if(knob) then
2315	localmaster=master
2316	call ass(dtant)
2317	call varfk1(S1)
2318	dtant%t= tan(varf1)
2319	master=localmaster
2320	else
2321	dtant%r= tan(S1%r)
2322	dtant%kind=1
2323	endif
2324	case default
2325	w_p=0
2326	w_p%nc=2
2327	w_p%fc='((1X,A72,/,1x,a72))'
2328	w_p%fi='(2((1X,i4)))'
2329	w_p%c(1)= " trouble in dtant "
2330	w_p%c(2)= "s1%kind "
2331	w_p=(/s1%kind /)
2332	! call !write_e(0)
2333	end select
2334	END FUNCTION dtant
2335
2336	FUNCTION dtandt( S1 )
2337	implicit none
2338	TYPE (real_8) dtandt
2339	TYPE (real_8), INTENT (IN) :: S1
2340	integer localmaster
2341
2342	select case(s1%kind)
2343	case(m1)
2344	dtandt%r=tan(s1%r*DEG_TO_RAD_)
2345	dtandt%kind=1
2346	case(m2)
2347	localmaster=master
2348	call ass(dtandt)
2349	dtandt%t= tan(s1%t*DEG_TO_RAD_)
2350	master=localmaster
2351	case(m3)
2352	if(knob) then
2353	localmaster=master
2354	call ass(dtandt)
2355	call varfk1(S1)
2356	dtandt%t= tan(varf1*DEG_TO_RAD_)
2357	master=localmaster
2358	else
2359	dtandt%r= tan(S1%r*DEG_TO_RAD_)
2360	dtandt%kind=1
2361	endif
2362	case default
2363	w_p=0
2364	w_p%nc=2
2365	w_p%fc='((1X,A72,/,1x,a72))'
2366	w_p%fi='(2((1X,i4)))'
2367	w_p%c(1)= " trouble in dtandt "
2368	w_p%c(2)= "s1%kind "
2369	w_p=(/s1%kind /)
2370	! call !write_e(0)
2371	end select
2372	END FUNCTION dtandt
2373
2374	FUNCTION absoftdtandr( S1 )
2375	implicit none
2376	real(dp) absoftdtandr
2377	real(dp), INTENT (IN) :: S1
2378	absoftdtandr=tan(s1*DEG_TO_RAD_)
2379	END FUNCTION absoftdtandr
2380
2381	FUNCTION dcost( S1 )
2382	implicit none
2383	TYPE (real_8) dcost
2384	TYPE (real_8), INTENT (IN) :: S1
2385	integer localmaster
2386
2387	select case(s1%kind)
2388	case(m1)
2389	dcost%r=cos(s1%r)
2390	dcost%kind=1
2391	case(m2)
2392	localmaster=master
2393	call ass(dcost)
2394	dcost%t= cos(s1%t)
2395	master=localmaster
2396	case(m3)
2397	if(knob) then
2398	localmaster=master
2399	call ass(dcost)
2400	call varfk1(S1)
2401	dcost%t= cos(varf1)
2402	master=localmaster
2403	else
2404	dcost%r= cos(S1%r)
2405	dcost%kind=1
2406	endif
2407	case default
2408	w_p=0
2409	w_p%nc=2
2410	w_p%fc='((1X,A72,/,1x,a72))'
2411	w_p%fi='(2((1X,i4)))'
2412	w_p%c(1)= " trouble in dcost "
2413	w_p%c(2)= "s1%kind "
2414	w_p=(/s1%kind /)
2415	! call !write_e(0)
2416	end select
2417	END FUNCTION dcost
2418
2419	FUNCTION dcosdt( S1 )
2420	implicit none
2421	TYPE (real_8) dcosdt
2422	TYPE (real_8), INTENT (IN) :: S1
2423	integer localmaster
2424
2425	select case(s1%kind)
2426	case(m1)
2427	dcosdt%r=cos(s1%r*DEG_TO_RAD_)
2428	dcosdt%kind=1
2429	case(m2)
2430	localmaster=master
2431	call ass(dcosdt)
2432	dcosdt%t= cos(s1%t*DEG_TO_RAD_)
2433	master=localmaster
2434	case(m3)
2435	if(knob) then
2436	localmaster=master
2437	call ass(dcosdt)
2438	call varfk1(S1)
2439	dcosdt%t= cos(varf1*DEG_TO_RAD_)
2440	master=localmaster
2441	else
2442	dcosdt%r= cos(S1%r)
2443	dcosdt%kind=1
2444	endif
2445	case default
2446	w_p=0
2447	w_p%nc=2
2448	w_p%fc='((1X,A72,/,1x,a72))'
2449	w_p%fi='(2((1X,i4)))'
2450	w_p%c(1)= " trouble in dcosdt "
2451	w_p%c(2)= "s1%kind "
2452	w_p=(/s1%kind /)
2453	! call !write_e(0)
2454	end select
2455	END FUNCTION dcosdt
2456
2457	FUNCTION absoftdcosdr( S1 )
2458	implicit none
2459	real(dp) absoftdcosdr
2460	real(dp) , INTENT (IN) :: S1
2461	absoftdcosdr=cos(s1*DEG_TO_RAD_ )
2462	END FUNCTION absoftdcosdr
2463
2464
2465	FUNCTION dsint( S1 )
2466	implicit none
2467	TYPE (real_8) dsint
2468	TYPE (real_8), INTENT (IN) :: S1
2469	integer localmaster
2470
2471	select case(s1%kind)
2472	case(m1)
2473	dsint%r=sin(s1%r)
2474	dsint%kind=1
2475	case(m2)
2476	localmaster=master
2477	call ass(dsint)
2478	dsint%t= sin(s1%t)
2479	master=localmaster
2480	case(m3)
2481	if(knob) then
2482	localmaster=master
2483	call ass(dsint)
2484	call varfk1(S1)
2485	dsint%t= sin(varf1)
2486	master=localmaster
2487	else
2488	dsint%r= sin(S1%r)
2489	dsint%kind=1
2490	endif
2491	case default
2492	w_p=0
2493	w_p%nc=2
2494	w_p%fc='((1X,A72,/,1x,a72))'
2495	w_p%fi='(2((1X,i4)))'
2496	w_p%c(1)= " trouble in dsint "
2497	w_p%c(2)= "s1%kind "
2498	w_p=(/s1%kind /)
2499	! call !write_e(0)
2500	end select
2501	END FUNCTION dsint
2502
2503	FUNCTION dsindt( S1 )
2504	implicit none
2505	TYPE (real_8) dsindt
2506	TYPE (real_8), INTENT (IN) :: S1
2507	integer localmaster
2508
2509	select case(s1%kind)
2510	case(m1)
2511	dsindt%r=sin(s1%r*DEG_TO_RAD_ )
2512	dsindt%kind=1
2513	case(m2)
2514	localmaster=master
2515	call ass(dsindt)
2516	dsindt%t= sin(s1%t*DEG_TO_RAD_)
2517	master=localmaster
2518	case(m3)
2519	if(knob) then
2520	localmaster=master
2521	call ass(dsindt)
2522	call varfk1(S1)
2523	dsindt%t= sin(varf1*DEG_TO_RAD_)
2524	master=localmaster
2525	else
2526	dsindt%r= sin(S1%r*DEG_TO_RAD_)
2527	dsindt%kind=1
2528	endif
2529	case default
2530	w_p=0
2531	w_p%nc=2
2532	w_p%fc='((1X,A72,/,1x,a72))'
2533	w_p%fi='(2((1X,i4)))'
2534	w_p%c(1)= " trouble in dsindt "
2535	w_p%c(2)= "s1%kind "
2536	w_p=(/s1%kind /)
2537	! call !write_e(0)
2538	end select
2539	END FUNCTION dsindt
2540
2541	FUNCTION absoftdsindr( S1 )
2542	implicit none
2543	real(dp) absoftdsindr
2544	real(dp) , INTENT (IN) :: S1
2545	absoftdsindr=sin(s1*DEG_TO_RAD_ )
2546	END FUNCTION absoftdsindr
2547
2548
2549	FUNCTION dlogt( S1 )
2550	implicit none
2551	TYPE (real_8) dlogt
2552	TYPE (real_8), INTENT (IN) :: S1
2553	integer localmaster
2554
2555	select case(s1%kind)
2556	case(m1)
2557	dlogt%r=loge(s1%r)
2558	dlogt%kind=1
2559	case(m2)
2560	localmaster=master
2561	call ass(dlogt)
2562	dlogt%t= log(s1%t)
2563	master=localmaster
2564	case(m3)
2565	if(knob) then
2566	localmaster=master
2567	call ass(dlogt)
2568	call varfk1(S1)
2569	dlogt%t= log(varf1)
2570	master=localmaster
2571	else
2572	dlogt%r= log(S1%r)
2573	dlogt%kind=1
2574	endif
2575	case default
2576	w_p=0
2577	w_p%nc=2
2578	w_p%fc='((1X,A72,/,1x,a72))'
2579	w_p%fi='(2((1X,i4)))'
2580	w_p%c(1)= " trouble in dlogt "
2581	w_p%c(2)= "s1%kind "
2582	w_p=(/s1%kind /)
2583	! call !write_e(0)
2584	end select
2585	END FUNCTION dlogt
2586
2587	FUNCTION dsqrtt( S1 )
2588	implicit none
2589	TYPE (real_8) dsqrtt
2590	TYPE (real_8), INTENT (IN) :: S1
2591	integer localmaster
2592
2593	select case(s1%kind)
2594	case(m1)
2595	! dsqrtt%r=SQRT(s1%r)
2596	dsqrtt%r=root(s1%r)
2597	dsqrtt%kind=1
2598	case(m2)
2599	localmaster=master
2600	call ass(dsqrtt)
2601	dsqrtt%t= SQRT(s1%t)
2602	master=localmaster
2603	case(m3)
2604	if(knob) then
2605	localmaster=master
2606	call ass(dsqrtt)
2607	call varfk1(S1)
2608	dsqrtt%t= SQRT(varf1)
2609	master=localmaster
2610	else
2611	! dsqrtt%r= SQRT(S1%r)
2612	dsqrtt%r= root(S1%r)
2613	dsqrtt%kind=1
2614	endif
2615	case default
2616	w_p=0
2617	w_p%nc=2
2618	w_p%fc='((1X,A72,/,1x,a72))'
2619	w_p%fi='(2((1X,i4)))'
2620	w_p%c(1)= " trouble in dsqrtt "
2621	w_p%c(2)= "s1%kind "
2622	w_p=(/s1%kind /)
2623	! call !write_e(0)
2624	end select
2625	END FUNCTION dsqrtt
2626
2627
2628
2629	FUNCTION POW( S1, S2 )
2630	implicit none
2631	TYPE (real_8) POW
2632	TYPE (real_8), INTENT (IN) :: S1
2633	integer , INTENT (IN) :: S2
2634	integer localmaster
2635
2636	select case(s1%kind)
2637	case(m1)
2638	POW%r=s1%r**s2
2639	POW%kind=1
2640	case(m2)
2641	localmaster=master
2642	call ass(POW)
2643	POW%t= s1%t**s2
2644	master=localmaster
2645	case(m3)
2646	if(knob) then
2647	localmaster=master
2648	call ass(POW)
2649	call varfk1(S1)
2650	POW%t= varf1**s2
2651	master=localmaster
2652	else
2653	POW%r= S1%r**s2
2654	POW%kind=1
2655	endif
2656	case default
2657	w_p=0
2658	w_p%nc=2
2659	w_p%fc='((1X,A72,/,1x,a72))'
2660	w_p%fi='(2((1X,i4)))'
2661	w_p%c(1)= " trouble in POW "
2662	w_p%c(2)= "s1%kind "
2663	w_p=(/s1%kind /)
2664	! call !write_e(0)
2665	end select
2666	END FUNCTION POW
2667
2668	FUNCTION POWR( S1, S2 )
2669	implicit none
2670	TYPE (real_8) POWR
2671	TYPE (real_8), INTENT (IN) :: S1
2672	REAL(SP) , INTENT (IN) :: S2
2673	integer localmaster
2674
2675	if(real_warning) call real_stop
2676	select case(s1%kind)
2677	case(m1)
2678	POWR%r=s1%r**REAL(s2,kind=DP)
2679	POWR%kind=1
2680	case(m2)
2681	localmaster=master
2682	call ass(POWR)
2683	POWR%t= s1%t**REAL(s2,kind=DP)
2684	master=localmaster
2685	case(m3)
2686	if(knob) then
2687	localmaster=master
2688	call ass(POWR)
2689	call varfk1(S1)
2690	POWR%t= varf1**REAL(s2,kind=DP)
2691	master=localmaster
2692	else
2693	POWR%r= S1%r**REAL(s2,kind=DP)
2694	POWR%kind=1
2695	endif
2696	case default
2697	w_p=0
2698	w_p%nc=2
2699	w_p%fc='((1X,A72,/,1x,a72))'
2700	w_p%fi='(2((1X,i4)))'
2701	w_p%c(1)= " trouble in POWR "
2702	w_p%c(2)= "s1%kind "
2703	w_p=(/s1%kind /)
2704	! call !write_e(0)
2705	end select
2706	END FUNCTION POWR
2707
2708	FUNCTION POWR8( S1, S2 )
2709	implicit none
2710	TYPE (real_8) POWR8
2711	TYPE (real_8), INTENT (IN) :: S1
2712	real(dp) , INTENT (IN) :: S2
2713	integer localmaster
2714
2715	select case(s1%kind)
2716	case(m1)
2717	POWR8%r=s1%r**s2
2718	POWR8%kind=1
2719	case(m2)
2720	localmaster=master
2721	call ass(POWR8)
2722	POWR8%t= s1%t**s2
2723	master=localmaster
2724	case(m3)
2725	if(knob) then
2726	localmaster=master
2727	call ass(POWR8)
2728	call varfk1(S1)
2729	POWR8%t= varf1**s2
2730	master=localmaster
2731	else
2732	POWR8%r= S1%r**s2
2733	POWR8%kind=1
2734	endif
2735	case default
2736	w_p=0
2737	w_p%nc=2
2738	w_p%fc='((1X,A72,/,1x,a72))'
2739	w_p%fi='(2((1X,i4)))'
2740	w_p%c(1)= " trouble in POWR8 "
2741	w_p%c(2)= "s1%kind "
2742	w_p=(/s1%kind /)
2743	! call !write_e(0)
2744	end select
2745	END FUNCTION POWR8
2746
2747
2748	FUNCTION unaryADD( S1 )
2749	implicit none
2750	TYPE (real_8) unaryADD
2751	TYPE (real_8), INTENT (IN) :: S1
2752	integer localmaster
2753
2754	select case(s1%kind)
2755	case(m1)
2756	unaryADD%r=s1%r
2757	unaryADD%kind=1
2758	case(m2)
2759	localmaster=master
2760	call ass(unaryADD)
2761	unaryADD%t= s1%t
2762	master=localmaster
2763	case(m3)
2764	if(knob) then
2765	localmaster=master
2766	call ass(unaryADD)
2767	call varfk1(S1)
2768	unaryADD%t= varf1
2769	master=localmaster
2770	else
2771	unaryADD%r= S1%r
2772	unaryADD%kind=1
2773	endif
2774	case default
2775	w_p=0
2776	w_p%nc=2
2777	w_p%fc='((1X,A72,/,1x,a72))'
2778	w_p%fi='(2((1X,i4)))'
2779	w_p%c(1)= " trouble in POWR8 "
2780	w_p%c(2)= "s1%kind "
2781	w_p=(/s1%kind /)
2782	! call !write_e(0)
2783	end select
2784	END FUNCTION unaryADD
2785
2786	FUNCTION unarySUB( S1 )
2787	implicit none
2788	TYPE (real_8) unarySUB
2789	TYPE (real_8), INTENT (IN) :: S1
2790	integer localmaster
2791
2792	select case(s1%kind)
2793	case(m1)
2794	unarySUB%r=-s1%r
2795	unarySUB%kind=1
2796	case(m2)
2797	localmaster=master
2798	call ass(unarySUB)
2799	unarySUB%t= -s1%t
2800	master=localmaster
2801	case(m3)
2802	if(knob) then
2803	localmaster=master
2804	call ass(unarySUB)
2805	call varfk1(S1)
2806	unarySUB%t= -varf1
2807	master=localmaster
2808	else
2809	unarySUB%r= -S1%r
2810	unarySUB%kind=1
2811	endif
2812	case default
2813	w_p=0
2814	w_p%nc=2
2815	w_p%fc='((1X,A72,/,1x,a72))'
2816	w_p%fi='(2((1X,i4)))'
2817	w_p%c(1)= " trouble in unarySUB "
2818	w_p%c(2)= "s1%kind "
2819	w_p=(/s1%kind /)
2820	! call !write_e(0)
2821	end select
2822	END FUNCTION unarySUB
2823
2824	FUNCTION add( S1, S2 )
2825	implicit none
2826	TYPE (real_8) add
2827	TYPE (real_8), INTENT (IN) :: S1, S2
2828	integer localmaster
2829
2830	select case(s1%kind+ms*s2%kind)
2831	case(m11)
2832	add%r=s1%r+s2%r
2833	add%kind=1
2834	case(m12,m21,m22)
2835	localmaster=master
2836	call ass(add)
2837	select case(s1%kind+ms*s2%kind)
2838	case(m21)
2839	add%t= s1%t+s2%r
2840	case(m12)
2841	add%t= s1%r+s2%t
2842	case(m22)
2843	add%t= s1%t+s2%t
2844	end select
2845	master=localmaster
2846
2847	case(m13,m31,m32,m23,m33)
2848	select case(s1%kind+ms*s2%kind)
2849	case(m31)
2850	if(knob) then
2851	localmaster=master
2852	call ass(add)
2853	call varfk1(S1)
2854	add%t= varf1+s2%r
2855	master=localmaster
2856	else
2857	add%r= s1%r+s2%r
2858	add%kind=1
2859	endif
2860	case(m13)
2861	if(knob) then
2862	localmaster=master
2863	call ass(add)
2864	call varfk1(S2)
2865	add%t= s1%r+varf1
2866	master=localmaster
2867	else
2868	add%r= s1%r+s2%r
2869	add%kind=1
2870	endif
2871	case(m32)
2872	localmaster=master
2873	call ass(add)
2874	if(knob) then
2875	call varfk1(S1)
2876	add%t= varf1+s2%t
2877	else
2878	add%t= s1%r+s2%t
2879	endif
2880	master=localmaster
2881	case(m23)
2882	localmaster=master
2883	call ass(add)
2884	if(knob) then
2885	call varfk1(S2)
2886	add%t= s1%t+varf1
2887	else
2888	add%t= s1%t+s2%r
2889	endif
2890	master=localmaster
2891	case(m33)
2892	if(knob) then
2893	localmaster=master
2894	call ass(add)
2895	call varfk1(S1)
2896	call varfk2(S2)
2897	add%t= varf1+varf2
2898	master=localmaster
2899	else
2900	add%r= s1%r+s2%r
2901	add%kind=1
2902	endif
2903	end select
2904	case default
2905	w_p=0
2906	w_p%nc=2
2907	w_p%fc='((1X,A72,/,1x,a72))'
2908	w_p%fi='(2((1X,i4)))'
2909	w_p%c(1)= " trouble in add "
2910	w_p%c(2)= "s1%kind ,s2%kind "
2911	w_p=(/s1%kind ,s2%kind/)
2912	! call !write_e(0)
2913	end select
2914	END FUNCTION add
2915
2916	FUNCTION subs( S1, S2 )
2917	implicit none
2918	TYPE (real_8) subs
2919	TYPE (real_8), INTENT (IN) :: S1, S2
2920	integer localmaster
2921
2922	select case(s1%kind+ms*s2%kind)
2923	case(m11)
2924	subs%r=s1%r-s2%r
2925	subs%kind=1
2926	case(m12,m21,m22)
2927	localmaster=master
2928	call ass(subs)
2929	select case(s1%kind+ms*s2%kind)
2930	case(m21)
2931	subs%t= s1%t-s2%r
2932	case(m12)
2933	subs%t= s1%r-s2%t
2934	case(m22)
2935	subs%t= s1%t-s2%t
2936	end select
2937	master=localmaster
2938	case(m13,m31,m32,m23,m33)
2939	select case(s1%kind+ms*s2%kind)
2940	case(m31)
2941
2942	if(knob) then
2943	localmaster=master
2944	call ass(subs)
2945	call varfk1(S1)
2946	subs%t= varf1-s2%r
2947	master=localmaster
2948	else
2949	subs%r= s1%r-s2%r
2950	subs%kind=1
2951	endif
2952	case(m13)
2953	if(knob) then
2954	localmaster=master
2955	call ass(subs)
2956	call varfk1(S2)
2957	subs%t= s1%r-varf1
2958	master=localmaster
2959	else
2960	subs%r= s1%r-s2%r
2961	subs%kind=1
2962	endif
2963	case(m32)
2964	localmaster=master
2965	call ass(subs)
2966	if(knob) then
2967	call varfk1(S1)
2968	subs%t= varf1-s2%t
2969	else
2970	subs%t= s1%r-s2%t
2971	endif
2972	master=localmaster
2973	case(m23)
2974	localmaster=master
2975	call ass(subs)
2976	if(knob) then
2977	call varfk1(S2)
2978	subs%t= s1%t-varf1
2979	else
2980	subs%t= s1%t-s2%r
2981	endif
2982	master=localmaster
2983	case(m33)
2984	if(knob) then
2985	localmaster=master
2986	call ass(subs)
2987	call varfk1(S1)
2988	call varfk2(S2)
2989	subs%t= varf1-varf2
2990	master=localmaster
2991	else
2992	subs%r= s1%r-s2%r
2993	subs%kind=1
2994	endif
2995	end select
2996	case default
2997	w_p=0
2998	w_p%nc=2
2999	w_p%fc='((1X,A72,/,1x,a72))'
3000	w_p%fi='(2((1X,i4)))'
3001	w_p%c(1)= " trouble in subs "
3002	w_p%c(2)= "s1%kind ,s2%kind "
3003	w_p=(/s1%kind ,s2%kind/)
3004	! call !write_e(0)
3005	end select
3006	END FUNCTION subs
3007
3008	FUNCTION daddsc( S1, S2 )
3009	implicit none
3010	TYPE (real_8) daddsc
3011	TYPE (real_8), INTENT (IN) :: S1
3012	real(dp) , INTENT (IN) :: S2
3013	integer localmaster
3014
3015	select case(s1%kind)
3016	case(m1)
3017	daddsc%r=s1%r+s2
3018	daddsc%kind=1
3019	case(m2)
3020	localmaster=master
3021	call ass(daddsc)
3022	daddsc%t= s1%t+s2
3023	master=localmaster
3024	case(m3)
3025	if(knob) then
3026	localmaster=master
3027	call ass(daddsc)
3028	call varfk1(S1)
3029	daddsc%t= varf1+s2
3030	master=localmaster
3031	else
3032	daddsc%r= s1%r+s2
3033	daddsc%kind=1
3034	endif
3035	case default
3036	w_p=0
3037	w_p%nc=2
3038	w_p%fc='((1X,A72,/,1x,a72))'
3039	w_p%fi='(2((1X,i4)))'
3040	w_p%c(1)= " trouble in daddsc "
3041	w_p%c(2)= "s1%kind "
3042	w_p=(/s1%kind /)
3043	! call !write_e(0)
3044	end select
3045	END FUNCTION daddsc
3046
3047	FUNCTION dscadd( S2, S1 )
3048	implicit none
3049	TYPE (real_8) dscadd
3050	TYPE (real_8), INTENT (IN) :: S1
3051	real(dp) , INTENT (IN) :: S2
3052	integer localmaster
3053
3054	select case(s1%kind)
3055	case(m1)
3056	dscadd%r=s1%r+s2
3057	dscadd%kind=1
3058	case(m2)
3059	localmaster=master
3060	call ass(dscadd)
3061	dscadd%t= s1%t+s2
3062	master=localmaster
3063	case(m3)
3064	if(knob) then
3065	localmaster=master
3066	call ass(dscadd)
3067	call varfk1(S1)
3068	dscadd%t= varf1+s2
3069	master=localmaster
3070	else
3071	dscadd%r= S1%r+s2
3072	dscadd%kind=1
3073	endif
3074	case default
3075	w_p=0
3076	w_p%nc=2
3077	w_p%fc='((1X,A72,/,1x,a72))'
3078	w_p%fi='(2((1X,i4)))'
3079	w_p%c(1)= " trouble in dscadd "
3080	w_p%c(2)= "s1%kind "
3081	w_p=(/s1%kind /)
3082	! call !write_e(0)
3083	end select
3084	END FUNCTION dscadd
3085
3086	FUNCTION dsubsc( S1, S2 )
3087	implicit none
3088	TYPE (real_8) dsubsc
3089	TYPE (real_8), INTENT (IN) :: S1
3090	real(dp) , INTENT (IN) :: S2
3091	integer localmaster
3092
3093	select case(s1%kind)
3094	case(m1)
3095	dsubsc%r=s1%r-s2
3096	dsubsc%kind=1
3097	case(m2)
3098	localmaster=master
3099	call ass(dsubsc)
3100	dsubsc%t= s1%t-s2
3101	master=localmaster
3102	case(m3)
3103	if(knob) then
3104	localmaster=master
3105	call ass(dsubsc)
3106	call varfk1(S1)
3107	dsubsc%t= varf1-s2
3108	master=localmaster
3109	else
3110	dsubsc%r=s1%r-s2
3111	dsubsc%kind=1
3112	endif
3113	case default
3114	w_p=0
3115	w_p%nc=2
3116	w_p%fc='((1X,A72,/,1x,a72))'
3117	w_p%fi='(2((1X,i4)))'
3118	w_p%c(1)= " trouble in dsubsc "
3119	w_p%c(2)= "s1%kind "
3120	w_p=(/s1%kind /)
3121	! call !write_e(0)
3122	end select
3123	END FUNCTION dsubsc
3124
3125	FUNCTION dscsub( S2, S1 )
3126	implicit none
3127	TYPE (real_8) dscsub
3128	TYPE (real_8), INTENT (IN) :: S1
3129	real(dp) , INTENT (IN) :: S2
3130	integer localmaster
3131
3132	select case(s1%kind)
3133	case(m1)
3134	dscsub%r=s2-s1%r
3135	dscsub%kind=1
3136	case(m2)
3137	localmaster=master
3138	call ass(dscsub)
3139	dscsub%t=s2-s1%t
3140	master=localmaster
3141	case(m3)
3142	if(knob) then
3143	localmaster=master
3144	call ass(dscsub)
3145	call varfk1(S1)
3146	dscsub%t=s2-varf1
3147	master=localmaster
3148	else
3149	dscsub%r=s2-s1%r
3150	dscsub%kind=1
3151	endif
3152	case default
3153	w_p=0
3154	w_p%nc=2
3155	w_p%fc='((1X,A72,/,1x,a72))'
3156	w_p%fi='(2((1X,i4)))'
3157	w_p%c(1)= " trouble in dscsub "
3158	w_p%c(2)= "s1%kind "
3159	w_p=(/s1%kind /)
3160	! call !write_e(0)
3161	end select
3162	END FUNCTION dscsub
3163
3164	FUNCTION addsc( S1, S2 )
3165	implicit none
3166	TYPE (real_8) addsc
3167	TYPE (real_8), INTENT (IN) :: S1
3168	real(sp) , INTENT (IN) :: S2
3169	integer localmaster
3170
3171	if(real_warning) call real_stop
3172	select case(s1%kind)
3173	case(m1)
3174	addsc%r=s1%r+REAL(s2,kind=DP)
3175	addsc%kind=1
3176	case(m2)
3177	localmaster=master
3178	call ass(addsc)
3179	addsc%t= s1%t+REAL(s2,kind=DP)
3180	master=localmaster
3181	case(m3)
3182	if(knob) then
3183	localmaster=master
3184	call ass(addsc)
3185	call varfk1(S1)
3186	addsc%t= varf1+REAL(s2,kind=DP)
3187	master=localmaster
3188	else
3189	addsc%r=s1%r+REAL(s2,kind=DP)
3190	addsc%kind=1
3191	endif
3192	case default
3193	w_p=0
3194	w_p%nc=2
3195	w_p%fc='((1X,A72,/,1x,a72))'
3196	w_p%fi='(2((1X,i4)))'
3197	w_p%c(1)= " trouble in addsc "
3198	w_p%c(2)= "s1%kind "
3199	w_p=(/s1%kind /)
3200	! call !write_e(0)
3201	end select
3202	END FUNCTION addsc
3203
3204	FUNCTION scadd( S2, S1 )
3205	implicit none
3206	TYPE (real_8) scadd
3207	TYPE (real_8), INTENT (IN) :: S1
3208	real(sp), INTENT (IN) :: S2
3209	integer localmaster
3210
3211	if(real_warning) call real_stop
3212	select case(s1%kind)
3213	case(m1)
3214	scadd%r=s1%r+REAL(s2,kind=DP)
3215	scadd%kind=1
3216	case(m2)
3217	localmaster=master
3218	call ass(scadd)
3219	scadd%t= s1%t+REAL(s2,kind=DP)
3220	master=localmaster
3221	case(m3)
3222	if(knob) then
3223	localmaster=master
3224	call ass(scadd)
3225	call varfk1(S1)
3226	scadd%t= varf1+REAL(s2,kind=DP)
3227	master=localmaster
3228	else
3229	scadd%r=s1%r+REAL(s2,kind=DP)
3230	scadd%kind=1
3231	endif
3232	case default
3233	w_p=0
3234	w_p%nc=2
3235	w_p%fc='((1X,A72,/,1x,a72))'
3236	w_p%fi='(2((1X,i4)))'
3237	w_p%c(1)= " trouble in addsc "
3238	w_p%c(2)= "s1%kind "
3239	w_p=(/s1%kind /)
3240	! call !write_e(0)
3241	end select
3242	END FUNCTION scadd
3243
3244	FUNCTION subsc( S1, S2 )
3245	implicit none
3246	TYPE (real_8) subsc
3247	TYPE (real_8), INTENT (IN) :: S1
3248	real(sp) , INTENT (IN) :: S2
3249	integer localmaster
3250
3251	if(real_warning) call real_stop
3252	select case(s1%kind)
3253	case(m1)
3254	subsc%r=s1%r-REAL(s2,kind=DP)
3255	subsc%kind=1
3256	case(m2)
3257	localmaster=master
3258	call ass(subsc)
3259	subsc%t= s1%t-REAL(s2,kind=DP)
3260	master=localmaster
3261	case(m3)
3262	if(knob) then
3263	localmaster=master
3264	call ass(subsc)
3265	call varfk1(S1)
3266	subsc%t= varf1-REAL(s2,kind=DP)
3267	master=localmaster
3268	else
3269	subsc%r=s1%r-REAL(s2,kind=DP)
3270	subsc%kind=1
3271	endif
3272	case default
3273	w_p=0
3274	w_p%nc=2
3275	w_p%fc='((1X,A72,/,1x,a72))'
3276	w_p%fi='(2((1X,i4)))'
3277	w_p%c(1)= " trouble in subsc "
3278	w_p%c(2)= "s1%kind "
3279	w_p=(/s1%kind /)
3280	! call !write_e(0)
3281	end select
3282	END FUNCTION subsc
3283
3284	FUNCTION scsub( S2, S1 )
3285	implicit none
3286	TYPE (real_8) scsub
3287	TYPE (real_8), INTENT (IN) :: S1
3288	real(sp) , INTENT (IN) :: S2
3289	integer localmaster
3290
3291	if(real_warning) call real_stop
3292	select case(s1%kind)
3293	case(m1)
3294	scsub%r=REAL(s2,kind=DP)-s1%r
3295	scsub%kind=1
3296	case(m2)
3297	localmaster=master
3298	call ass(scsub)
3299	scsub%t=REAL(s2,kind=DP)-s1%t
3300	master=localmaster
3301	case(m3)
3302	if(knob) then
3303	localmaster=master
3304	call ass(scsub)
3305	call varfk1(S1)
3306	scsub%t=REAL(s2,kind=DP)-varf1
3307	master=localmaster
3308	else
3309	scsub%r=REAL(s2,kind=DP)-s1%r
3310	scsub%kind=1
3311	endif
3312	case default
3313	w_p=0
3314	w_p%nc=2
3315	w_p%fc='((1X,A72,/,1x,a72))'
3316	w_p%fi='(2((1X,i4)))'
3317	w_p%c(1)= " trouble in scsub "
3318	w_p%c(2)= "s1%kind "
3319	w_p=(/s1%kind /)
3320	! call !write_e(0)
3321	end select
3322	END FUNCTION scsub
3323
3324	FUNCTION iaddsc( S1, S2 )
3325	implicit none
3326	TYPE (real_8) iaddsc
3327	TYPE (real_8), INTENT (IN) :: S1
3328	integer , INTENT (IN) :: S2
3329	integer localmaster
3330
3331	select case(s1%kind)
3332	case(m1)
3333	iaddsc%r=s1%r+REAL(s2,kind=DP)
3334	iaddsc%kind=1
3335	case(m2)
3336	localmaster=master
3337	call ass(iaddsc)
3338	iaddsc%t= s1%t+REAL(s2,kind=DP)
3339	master=localmaster
3340	case(m3)
3341	if(knob) then
3342	localmaster=master
3343	call ass(iaddsc)
3344	call varfk1(S1)
3345	iaddsc%t= varf1+REAL(s2,kind=DP)
3346	master=localmaster
3347	else
3348	iaddsc%r=s1%r+REAL(s2,kind=DP)
3349	iaddsc%kind=1
3350	endif
3351	case default
3352	w_p=0
3353	w_p%nc=2
3354	w_p%fc='((1X,A72,/,1x,a72))'
3355	w_p%fi='(2((1X,i4)))'
3356	w_p%c(1)= " trouble in iaddsc "
3357	w_p%c(2)= "s1%kind "
3358	w_p=(/s1%kind /)
3359	! call !write_e(0)
3360	end select
3361	END FUNCTION iaddsc
3362
3363	FUNCTION iscadd( S2, S1 )
3364	implicit none
3365	TYPE (real_8) iscadd
3366	TYPE (real_8), INTENT (IN) :: S1
3367	integer, INTENT (IN) :: S2
3368	integer localmaster
3369
3370	select case(s1%kind)
3371	case(m1)
3372	iscadd%r=s1%r+REAL(s2,kind=DP)
3373	iscadd%kind=1
3374	case(m2)
3375	localmaster=master
3376	call ass(iscadd)
3377	iscadd%t= s1%t+REAL(s2,kind=DP)
3378	master=localmaster
3379	case(m3)
3380	if(knob) then
3381	localmaster=master
3382	call ass(iscadd)
3383	call varfk1(S1)
3384	iscadd%t= varf1+REAL(s2,kind=DP)
3385	master=localmaster
3386	else
3387	iscadd%r=s1%r+REAL(s2,kind=DP)
3388	iscadd%kind=1
3389	endif
3390	case default
3391	w_p=0
3392	w_p%nc=2
3393	w_p%fc='((1X,A72,/,1x,a72))'
3394	w_p%fi='(2((1X,i4)))'
3395	w_p%c(1)= " trouble in iscadd "
3396	w_p%c(2)= "s1%kind "
3397	w_p=(/s1%kind /)
3398	! call !write_e(0)
3399	end select
3400	END FUNCTION iscadd
3401
3402	FUNCTION isubsc( S1, S2 )
3403	implicit none
3404	TYPE (real_8) isubsc
3405	TYPE (real_8), INTENT (IN) :: S1
3406	integer , INTENT (IN) :: S2
3407	integer localmaster
3408
3409	select case(s1%kind)
3410	case(m1)
3411	isubsc%r=s1%r-REAL(s2,kind=DP)
3412	isubsc%kind=1
3413	case(m2)
3414	localmaster=master
3415	call ass(isubsc)
3416	isubsc%t= s1%t-REAL(s2,kind=DP)
3417	master=localmaster
3418	case(m3)
3419	if(knob) then
3420	localmaster=master
3421	call ass(isubsc)
3422	call varfk1(S1)
3423	isubsc%t= varf1-REAL(s2,kind=DP)
3424	master=localmaster
3425	else
3426	isubsc%r=s1%r-REAL(s2,kind=DP)
3427	isubsc%kind=1
3428	endif
3429	case default
3430	w_p=0
3431	w_p%nc=2
3432	w_p%fc='((1X,A72,/,1x,a72))'
3433	w_p%fi='(2((1X,i4)))'
3434	w_p%c(1)= " trouble in isubsc "
3435	w_p%c(2)= "s1%kind "
3436	w_p=(/s1%kind /)
3437	! call !write_e(0)
3438	end select
3439	END FUNCTION isubsc
3440
3441	FUNCTION iscsub( S2, S1 )
3442	implicit none
3443	TYPE (real_8) iscsub
3444	TYPE (real_8), INTENT (IN) :: S1
3445	integer , INTENT (IN) :: S2
3446	integer localmaster
3447
3448	select case(s1%kind)
3449	case(m1)
3450	iscsub%r=REAL(s2,kind=DP)-s1%r
3451	iscsub%kind=1
3452	case(m2)
3453	localmaster=master
3454	call ass(iscsub)
3455	iscsub%t=REAL(s2,kind=DP)-s1%t
3456	master=localmaster
3457	case(m3)
3458	if(knob) then
3459	localmaster=master
3460	call ass(iscsub)
3461	call varfk1(S1)
3462	iscsub%t=REAL(s2,kind=DP)-varf1
3463	master=localmaster
3464	else
3465	iscsub%r=REAL(s2,kind=DP)-s1%r
3466	iscsub%kind=1
3467	endif
3468	case default
3469	w_p=0
3470	w_p%nc=2
3471	w_p%fc='((1X,A72,/,1x,a72))'
3472	w_p%fi='(2((1X,i4)))'
3473	w_p%c(1)= " trouble in iscsub "
3474	w_p%c(2)= "s1%kind "
3475	w_p=(/s1%kind /)
3476	! call !write_e(0)
3477	end select
3478	END FUNCTION iscsub
3479
3480	FUNCTION mul( S1, S2 )
3481	implicit none
3482	TYPE (real_8) mul
3483	TYPE (real_8), INTENT (IN) :: S1, S2
3484	integer localmaster
3485
3486	select case(s1%kind+ms*s2%kind)
3487	case(m11)
3488	mul%r=s1%r*s2%r
3489	mul%kind=1
3490	case(m12,m21,m22)
3491	localmaster=master
3492	call ass(mul)
3493	select case(s1%kind+ms*s2%kind)
3494	case(m21)
3495	mul%t= s1%t*s2%r
3496	case(m12)
3497	mul%t= s1%r*s2%t
3498	case(m22)
3499	mul%t= s1%t*s2%t
3500	end select
3501	master=localmaster
3502
3503	case(m13,m31,m32,m23,m33)
3504	select case(s1%kind+ms*s2%kind)
3505	case(m31)
3506
3507	if(knob) then
3508	localmaster=master
3509	call ass(mul)
3510	call varfk1(S1)
3511	mul%t= varf1*s2%r
3512	master=localmaster
3513	else
3514	mul%r= s1%r*s2%r
3515	mul%kind=1
3516	endif
3517	case(m13)
3518	if(knob) then
3519	localmaster=master
3520	call ass(mul)
3521	call varfk1(S2)
3522	mul%t= s1%r*varf1
3523	master=localmaster
3524	else
3525	mul%r= s1%r*s2%r
3526	mul%kind=1
3527	endif
3528	case(m32)
3529	localmaster=master
3530	call ass(mul)
3531	if(knob) then
3532	call varfk1(S1)
3533	mul%t= varf1*s2%t
3534	else
3535	mul%t= s1%r*s2%t
3536	endif
3537	master=localmaster
3538	case(m23)
3539	localmaster=master
3540	call ass(mul)
3541	if(knob) then
3542	call varfk1(S2)
3543	mul%t= s1%t*varf1
3544	else
3545	mul%t= s1%t*s2%r
3546	endif
3547	master=localmaster
3548	case(m33)
3549	if(knob) then
3550	localmaster=master
3551	call ass(mul)
3552	call varfk1(S1)
3553	call varfk2(S2)
3554	mul%t= varf1*varf2
3555	master=localmaster
3556	else
3557	mul%r= s1%r*s2%r
3558	mul%kind=1
3559	endif
3560	end select
3561	case default
3562	w_p=0
3563	w_p%nc=2
3564	w_p%fc='((1X,A72,/,1x,a72))'
3565	w_p%fi='(2((1X,i4)))'
3566	w_p%c(1)= " trouble in mul "
3567	w_p%c(2)= "s1%kind ,s2%kind "
3568	w_p=(/s1%kind ,s2%kind/)
3569	! call !write_e(0)
3570	end select
3571	END FUNCTION mul
3572
3573	FUNCTION div( S1, S2 )
3574	implicit none
3575	TYPE (real_8) div
3576	TYPE (real_8), INTENT (IN) :: S1, S2
3577	integer localmaster
3578
3579	select case(s1%kind+ms*s2%kind)
3580	case(m11)
3581	div%r=s1%r/s2%r
3582	div%kind=1
3583	case(m12,m21,m22)
3584	localmaster=master
3585	call ass(div)
3586	select case(s1%kind+ms*s2%kind)
3587	case(m21)
3588	div%t= s1%t/s2%r
3589	case(m12)
3590	div%t= s1%r/s2%t
3591	case(m22)
3592	div%t= s1%t/s2%t
3593	end select
3594	master=localmaster
3595	case(m13,m31,m32,m23,m33)
3596	select case(s1%kind+ms*s2%kind)
3597	case(m31)
3598
3599	if(knob) then
3600	localmaster=master
3601	call ass(div)
3602	call varfk1(S1)
3603	div%t= varf1/s2%r
3604	master=localmaster
3605	else
3606	div%r= s1%r/s2%r
3607	div%kind=1
3608	endif
3609	case(m13)
3610	if(knob) then
3611	localmaster=master
3612	call ass(div)
3613	call varfk1(S2)
3614	div%t= s1%r/varf1
3615	master=localmaster
3616	else
3617	div%r= s1%r/s2%r
3618	div%kind=1
3619	endif
3620	case(m32)
3621	localmaster=master
3622	call ass(div)
3623	if(knob) then
3624	call varfk1(S1)
3625	div%t= varf1/s2%t
3626	else
3627	div%t= s1%r/s2%t
3628	endif
3629	master=localmaster
3630	case(m23)
3631	localmaster=master
3632	call ass(div)
3633	if(knob) then
3634	call varfk1(S2)
3635	div%t= s1%t/varf1
3636	else
3637	div%t= s1%t/s2%r
3638	endif
3639	master=localmaster
3640	case(m33)
3641	if(knob) then
3642	localmaster=master
3643	call ass(div)
3644	call varfk1(S1)
3645	call varfk2(S2)
3646	div%t= varf1/varf2
3647	master=localmaster
3648	else
3649	div%r= s1%r/s2%r
3650	div%kind=1
3651	endif
3652	end select
3653	case default
3654	w_p=0
3655	w_p%nc=2
3656	w_p%fc='((1X,A72,/,1x,a72))'
3657	w_p%fi='(2((1X,i4)))'
3658	w_p%c(1)= " trouble in div "
3659	w_p%c(2)= "s1%kind ,s2%kind "
3660	w_p=(/s1%kind ,s2%kind/)
3661	! call !write_e(0)
3662	end select
3663	END FUNCTION div
3664
3665	FUNCTION dmulmapconcat( S1, S2 )
3666	implicit none
3667	TYPE (real_8) dmulmapconcat
3668	TYPE (real_8), INTENT (IN) :: S1
3669	type(damap) , INTENT (IN) :: S2
3670	integer localmaster
3671
3672	select case(s1%kind)
3673	case(m1)
3674	dmulmapconcat%r=s1%r
3675	dmulmapconcat%kind=1
3676	case(m2)
3677	localmaster=master
3678	call ass(dmulmapconcat)
3679	dmulmapconcat%t= s1%t*s2
3680	master=localmaster
3681	case(m3)
3682	if(knob) then
3683	localmaster=master
3684	call ass(dmulmapconcat)
3685	call varfk1(S1)
3686	dmulmapconcat%t= varf1*s2
3687	master=localmaster
3688	else
3689	dmulmapconcat%r=s1%r
3690	dmulmapconcat%kind=1
3691	endif
3692	case default
3693	w_p=0
3694	w_p%nc=2
3695	w_p%fc='((1X,A72,/,1x,a72))'
3696	w_p%fi='(2((1X,i4)))'
3697	w_p%c(1)= " trouble in dmulmapconcat "
3698	w_p%c(2)= "s1%kind "
3699	w_p=(/s1%kind /)
3700	! call !write_e(0)
3701	end select
3702	END FUNCTION dmulmapconcat
3703
3704
3705	FUNCTION dmulsc( S1, S2 )
3706	implicit none
3707	TYPE (real_8) dmulsc
3708	TYPE (real_8), INTENT (IN) :: S1
3709	real(dp) , INTENT (IN) :: S2
3710	integer localmaster
3711
3712	select case(s1%kind)
3713	case(m1)
3714	dmulsc%r=s1%r*s2
3715	dmulsc%kind=1
3716	case(m2)
3717	localmaster=master
3718	call ass(dmulsc)
3719	dmulsc%t= s1%t*s2
3720	master=localmaster
3721	case(m3)
3722	if(knob) then
3723	localmaster=master
3724	call ass(dmulsc)
3725	call varfk1(S1)
3726	dmulsc%t= varf1*s2
3727	master=localmaster
3728	else
3729	dmulsc%r=s1%r*s2
3730	dmulsc%kind=1
3731	endif
3732	case default
3733	w_p=0
3734	w_p%nc=2
3735	w_p%fc='((1X,A72,/,1x,a72))'
3736	w_p%fi='(2((1X,i4)))'
3737	w_p%c(1)= " trouble in dmulsc "
3738	w_p%c(2)= "s1%kind "
3739	w_p=(/s1%kind /)
3740	! call !write_e(0)
3741	end select
3742	END FUNCTION dmulsc
3743
3744	FUNCTION dscmul( S2, S1 )
3745	implicit none
3746	TYPE (real_8) dscmul
3747	TYPE (real_8), INTENT (IN) :: S1
3748	real(dp) , INTENT (IN) :: S2
3749	integer localmaster
3750
3751	select case(s1%kind)
3752	case(m1)
3753	dscmul%r=s1%r*s2
3754	dscmul%kind=1
3755	case(m2)
3756	localmaster=master
3757	call ass(dscmul)
3758	dscmul%t= s1%t*s2
3759	master=localmaster
3760	case(m3)
3761	if(knob) then
3762	localmaster=master
3763	call ass(dscmul)
3764	call varfk1(S1)
3765	dscmul%t= varf1*s2
3766	master=localmaster
3767	else
3768	dscmul%r=s1%r*s2
3769	dscmul%kind=1
3770	endif
3771	end select
3772	END FUNCTION dscmul
3773
3774	FUNCTION ddivsc( S1, S2 )
3775	implicit none
3776	TYPE (real_8) ddivsc
3777	TYPE (real_8), INTENT (IN) :: S1
3778	real(dp) , INTENT (IN) :: S2
3779	integer localmaster
3780
3781	select case(s1%kind)
3782	case(m1)
3783	ddivsc%r=s1%r/s2
3784	ddivsc%kind=1
3785	case(m2)
3786	localmaster=master
3787	call ass(ddivsc)
3788	ddivsc%t= s1%t/s2
3789	master=localmaster
3790	case(m3)
3791	if(knob) then
3792	localmaster=master
3793	call ass(ddivsc)
3794	call varfk1(S1)
3795	ddivsc%t= varf1/s2
3796	master=localmaster
3797	else
3798	ddivsc%r=s1%r/s2
3799	ddivsc%kind=1
3800	endif
3801	case default
3802	w_p=0
3803	w_p%nc=2
3804	w_p%fc='((1X,A72,/,1x,a72))'
3805	w_p%fi='(2((1X,i4)))'
3806	w_p%c(1)= " trouble in ddivsc "
3807	w_p%c(2)= "s1%kind "
3808	w_p=(/s1%kind /)
3809	! call !write_e(0)
3810	end select
3811	END FUNCTION ddivsc
3812
3813	FUNCTION dscdiv( S2, S1 )
3814	implicit none
3815	TYPE (real_8) dscdiv
3816	TYPE (real_8), INTENT (IN) :: S1
3817	real(dp) , INTENT (IN) :: S2
3818	integer localmaster
3819
3820	select case(s1%kind)
3821	case(m1)
3822	dscdiv%r=s2/s1%r
3823	dscdiv%kind=1
3824	case(m2)
3825	localmaster=master
3826	call ass(dscdiv)
3827	dscdiv%t= s2/s1%t
3828	master=localmaster
3829	case(m3)
3830	if(knob) then
3831	localmaster=master
3832	call ass(dscdiv)
3833	call varfk1(S1)
3834	dscdiv%t= s2/varf1
3835	master=localmaster
3836	else
3837	dscdiv%r=s2/s1%r
3838	dscdiv%kind=1
3839	endif
3840	case default
3841	w_p=0
3842	w_p%nc=2
3843	w_p%fc='((1X,A72,/,1x,a72))'
3844	w_p%fi='(2((1X,i4)))'
3845	w_p%c(1)= " trouble in ddivsc "
3846	w_p%c(2)= "s1%kind "
3847	w_p=(/s1%kind /)
3848	! call !write_e(0)
3849	end select
3850	END FUNCTION dscdiv
3851
3852	FUNCTION mulsc( S1, S2 )
3853	implicit none
3854	TYPE (real_8) mulsc
3855	TYPE (real_8), INTENT (IN) :: S1
3856	real(sp) , INTENT (IN) :: S2
3857	integer localmaster
3858
3859	if(real_warning) call real_stop
3860	select case(s1%kind)
3861	case(m1)
3862	mulsc%r=s1%r*REAL(s2,kind=DP)
3863	mulsc%kind=1
3864	case(m2)
3865	localmaster=master
3866	call ass(mulsc)
3867	mulsc%t= s1%t*REAL(s2,kind=DP)
3868	master=localmaster
3869	case(m3)
3870	if(knob) then
3871	localmaster=master
3872	call ass(mulsc)
3873	call varfk1(S1)
3874	mulsc%t= varf1*REAL(s2,kind=DP)
3875	master=localmaster
3876	else
3877	mulsc%r=s1%r*REAL(s2,kind=DP)
3878	mulsc%kind=1
3879	endif
3880	case default
3881	w_p=0
3882	w_p%nc=2
3883	w_p%fc='((1X,A72,/,1x,a72))'
3884	w_p%fi='(2((1X,i4)))'
3885	w_p%c(1)= " trouble in mulsc "
3886	w_p%c(2)= "s1%kind "
3887	w_p=(/s1%kind /)
3888	! call !write_e(0)
3889	end select
3890	END FUNCTION mulsc
3891
3892	FUNCTION scmul( S2, S1 )
3893	implicit none
3894	TYPE (real_8) scmul
3895	TYPE (real_8), INTENT (IN) :: S1
3896	real(sp) , INTENT (IN) :: S2
3897	integer localmaster
3898
3899	if(real_warning) call real_stop
3900	select case(s1%kind)
3901	case(m1)
3902	scmul%r=s1%r*REAL(s2,kind=DP)
3903	scmul%kind=1
3904	case(m2)
3905	localmaster=master
3906	call ass(scmul)
3907	scmul%t= s1%t*REAL(s2,kind=DP)
3908	master=localmaster
3909	case(m3)
3910	if(knob) then
3911	localmaster=master
3912	call ass(scmul)
3913	call varfk1(S1)
3914	scmul%t= varf1*REAL(s2,kind=DP)
3915	master=localmaster
3916	else
3917	scmul%r=s1%r*REAL(s2,kind=DP)
3918	scmul%kind=1
3919	endif
3920	case default
3921	w_p=0
3922	w_p%nc=2
3923	w_p%fc='((1X,A72,/,1x,a72))'
3924	w_p%fi='(2((1X,i4)))'
3925	w_p%c(1)= " trouble in scmul "
3926	w_p%c(2)= "s1%kind "
3927	w_p=(/s1%kind /)
3928	! call !write_e(0)
3929	end select
3930	END FUNCTION scmul
3931
3932	FUNCTION divsc( S1, S2 )
3933	implicit none
3934	TYPE (real_8) divsc
3935	TYPE (real_8), INTENT (IN) :: S1
3936	real(sp) , INTENT (IN) :: S2
3937	integer localmaster
3938
3939	if(real_warning) call real_stop
3940	select case(s1%kind)
3941	case(m1)
3942	divsc%r=s1%r/REAL(s2,kind=DP)
3943	divsc%kind=1
3944	case(m2)
3945	localmaster=master
3946	call ass(divsc)
3947	divsc%t= s1%t/REAL(s2,kind=DP)
3948	master=localmaster
3949	case(m3)
3950	if(knob) then
3951	localmaster=master
3952	call ass(divsc)
3953	call varfk1(S1)
3954	divsc%t= varf1/REAL(s2,kind=DP)
3955	master=localmaster
3956	else
3957	divsc%r=s1%r/REAL(s2,kind=DP)
3958	divsc%kind=1
3959	endif
3960	case default
3961	w_p=0
3962	w_p%nc=2
3963	w_p%fc='((1X,A72,/,1x,a72))'
3964	w_p%fi='(2((1X,i4)))'
3965	w_p%c(1)= " trouble in divsc "
3966	w_p%c(2)= "s1%kind "
3967	w_p=(/s1%kind /)
3968	! call !write_e(0)
3969	end select
3970	END FUNCTION divsc
3971
3972	FUNCTION scdiv( S2, S1 )
3973	implicit none
3974	TYPE (real_8) scdiv
3975	TYPE (real_8), INTENT (IN) :: S1
3976	real(sp) , INTENT (IN) :: S2
3977	integer localmaster
3978
3979	if(real_warning) call real_stop
3980	select case(s1%kind)
3981	case(m1)
3982	scdiv%r=REAL(s2,kind=DP)/s1%r
3983	scdiv%kind=1
3984	case(m2)
3985	localmaster=master
3986	call ass(scdiv)
3987	scdiv%t= REAL(s2,kind=DP)/s1%t
3988	master=localmaster
3989	case(m3)
3990	if(knob) then
3991	localmaster=master
3992	call ass(scdiv)
3993	call varfk1(S1)
3994	scdiv%t= REAL(s2,kind=DP)/varf1
3995	master=localmaster
3996	else
3997	scdiv%r=REAL(s2,kind=DP)/s1%r
3998	scdiv%kind=1
3999	endif
4000	case default
4001	w_p=0
4002	w_p%nc=2
4003	w_p%fc='((1X,A72,/,1x,a72))'
4004	w_p%fi='(2((1X,i4)))'
4005	w_p%c(1)= " trouble in scdiv "
4006	w_p%c(2)= "s1%kind "
4007	w_p=(/s1%kind /)
4008	! call !write_e(0)
4009	end select
4010	END FUNCTION scdiv
4011
4012	FUNCTION imulsc( S1, S2 )
4013	implicit none
4014	TYPE (real_8) imulsc
4015	TYPE (real_8), INTENT (IN) :: S1
4016	integer , INTENT (IN) :: S2
4017	integer localmaster
4018
4019	select case(s1%kind)
4020	case(m1)
4021	imulsc%r=s1%r*REAL(s2,kind=DP)
4022	imulsc%kind=1
4023	case(m2)
4024	localmaster=master
4025	call ass(imulsc)
4026	imulsc%t= s1%t*REAL(s2,kind=DP)
4027	master=localmaster
4028	case(m3)
4029	if(knob) then
4030	localmaster=master
4031	call ass(imulsc)
4032	call varfk1(S1)
4033	imulsc%t= varf1*REAL(s2,kind=DP)
4034	master=localmaster
4035	else
4036	imulsc%r=s1%r*REAL(s2,kind=DP)
4037	imulsc%kind=1
4038	endif
4039	case default
4040	w_p=0
4041	w_p%nc=2
4042	w_p%fc='((1X,A72,/,1x,a72))'
4043	w_p%fi='(2((1X,i4)))'
4044	w_p%c(1)= " trouble in imulsc "
4045	w_p%c(2)= "s1%kind "
4046	w_p=(/s1%kind /)
4047	! call !write_e(0)
4048	end select
4049	END FUNCTION imulsc
4050
4051	FUNCTION iscmul( S2, S1 )
4052	implicit none
4053	TYPE (real_8) iscmul
4054	TYPE (real_8), INTENT (IN) :: S1
4055	integer , INTENT (IN) :: S2
4056	integer localmaster
4057
4058	select case(s1%kind)
4059	case(m1)
4060	iscmul%r=s1%r*REAL(s2,kind=DP)
4061	iscmul%kind=1
4062	case(m2)
4063	localmaster=master
4064	call ass(iscmul)
4065	iscmul%t= s1%t*REAL(s2,kind=DP)
4066	master=localmaster
4067	case(m3)
4068	if(knob) then
4069	localmaster=master
4070	call ass(iscmul)
4071	call varfk1(S1)
4072	iscmul%t= varf1*REAL(s2,kind=DP)
4073	master=localmaster
4074	else
4075	iscmul%r=s1%r*REAL(s2,kind=DP)
4076	iscmul%kind=1
4077	endif
4078	case default
4079	w_p=0
4080	w_p%nc=2
4081	w_p%fc='((1X,A72,/,1x,a72))'
4082	w_p%fi='(2((1X,i4)))'
4083	w_p%c(1)= " trouble in iscmul "
4084	w_p%c(2)= "s1%kind "
4085	w_p=(/s1%kind /)
4086	! call !write_e(0)
4087	end select
4088	END FUNCTION iscmul
4089
4090	FUNCTION idivsc( S1, S2 )
4091	implicit none
4092	TYPE (real_8) idivsc
4093	TYPE (real_8), INTENT (IN) :: S1
4094	integer , INTENT (IN) :: S2
4095	integer localmaster
4096
4097	select case(s1%kind)
4098	case(m1)
4099	idivsc%r=s1%r/REAL(s2,kind=DP)
4100	idivsc%kind=1
4101	case(m2)
4102	localmaster=master
4103	call ass(idivsc)
4104	idivsc%t= s1%t/REAL(s2,kind=DP)
4105	master=localmaster
4106	case(m3)
4107	if(knob) then
4108	localmaster=master
4109	call ass(idivsc)
4110	call varfk1(S1)
4111	idivsc%t= varf1/REAL(s2,kind=DP)
4112	master=localmaster
4113	else
4114	idivsc%r=s1%r/REAL(s2,kind=DP)
4115	idivsc%kind=1
4116	endif
4117	case default
4118	w_p=0
4119	w_p%nc=2
4120	w_p%fc='((1X,A72,/,1x,a72))'
4121	w_p%fi='(2((1X,i4)))'
4122	w_p%c(1)= " trouble in idivsc "
4123	w_p%c(2)= "s1%kind "
4124	w_p=(/s1%kind /)
4125	! call !write_e(0)
4126	end select
4127	END FUNCTION idivsc
4128
4129	FUNCTION iscdiv( S2, S1 )
4130	implicit none
4131	TYPE (real_8) iscdiv
4132	TYPE (real_8), INTENT (IN) :: S1
4133	integer , INTENT (IN) :: S2
4134	integer localmaster
4135
4136	select case(s1%kind)
4137	case(m1)
4138	iscdiv%r=REAL(s2,kind=DP)/s1%r
4139	iscdiv%kind=1
4140	case(m2)
4141	localmaster=master
4142	call ass(iscdiv)
4143	iscdiv%t= REAL(s2,kind=DP)/s1%t
4144	master=localmaster
4145	case(m3)
4146	if(knob) then
4147	localmaster=master
4148	call ass(iscdiv)
4149	call varfk1(S1)
4150	iscdiv%t= REAL(s2,kind=DP)/varf1
4151	master=localmaster
4152	else
4153	iscdiv%r=REAL(s2,kind=DP)/s1%r
4154	iscdiv%kind=1
4155	endif
4156	case default
4157	w_p=0
4158	w_p%nc=2
4159	w_p%fc='((1X,A72,/,1x,a72))'
4160	w_p%fi='(2((1X,i4)))'
4161	w_p%c(1)= " trouble in iscdiv "
4162	w_p%c(2)= "s1%kind "
4163	w_p=(/s1%kind /)
4164	! call !write_e(0)
4165	end select
4166	END FUNCTION iscdiv
4167
4168	SUBROUTINE printpoly(S2,i,prec)
4169	implicit none
4170	integer ipause, mypauses
4171	type (real_8),INTENT(INOUT)::S2
4172	integer i
4173	real(dp), optional :: prec
4174
4175	if(s2%kind/=0) then
4176
4177	select case (s2%kind)
4178	case(m1)
4179	write(i,*) s2%r
4180	case(m2)
4181	call pri(S2%t,i,prec)
4182	case(m3)
4183	if(s2%i>0) then
4184	write(i,) s2%r," +",s2%s," x_",int(s2%i)
4185	else
4186	write(i,*) s2%r
4187	endif
4188	if(s2%alloc) then
4189	write(line,'(a41)') " weird Taylor part should be deallocated "
4190	ipause=mypauses(0,line)
4191	endif
4192	end select
4193	else
4194
4195	line= "Warning not defined in Printpoly (real polymorph)"
4196	ipause=mypauses(1,line)
4197
4198	endif
4199
4200	END SUBROUTINE printpoly
4201
4202	SUBROUTINE printdouble(S2,i)
4203	implicit none
4204	real(dp),INTENT(INOUT)::S2
4205	integer i
4206
4207	write(i,*) s2
4208
4209	END SUBROUTINE printdouble
4210
4211	SUBROUTINE printsingle(S2,i)
4212	implicit none
4213	real(sp),INTENT(INOUT)::S2
4214	integer i
4215
4216	write(i,*) s2
4217
4218	END SUBROUTINE printsingle
4219
4220	SUBROUTINE resetpoly(S2)
4221	implicit none
4222	type (real_8),INTENT(INOUT)::S2
4223
4224	if(s2%alloc) call killtpsa(s2%t)
4225	s2%alloc=f
4226	s2%kind=0
4227	s2%r=0.0_dp
4228	!s2%s=one
4229	!s2%i=0
4230
4231	END SUBROUTINE resetpoly
4232
4233
4234	SUBROUTINE resetpolyn(S2,K)
4235	implicit none
4236	type (real_8),INTENT(INOUT),dimension(:)::S2
4237	INTEGER,optional,INTENT(IN)::k
4238	INTEGER J,i,N
4239
4240	if(present(k)) then
4241	I=LBOUND(S2,DIM=1)
4242	N=LBOUND(S2,DIM=1)+K-1
4243	else
4244	I=LBOUND(S2,DIM=1)
4245	N=UBOUND(S2,DIM=1)
4246	endif
4247
4248	DO J=I,N
4249	call resetpoly(s2(j))
4250	enddo
4251
4252
4253	END SUBROUTINE resetpolyn
4254
4255
4256	SUBROUTINE resetpoly_R(S2,FL) ! STAYS REAL
4257	implicit none
4258	type (real_8),INTENT(INOUT)::S2
4259	logical(lp),INTENT(IN)::FL
4260
4261	if(s2%alloc) call killtpsa(s2%t)
4262	s2%alloc=f
4263	s2%kind=1
4264	s2%r=0.0_dp
4265
4266	IF(.NOT.FL) THEN
4267	s2%i=0
4268	s2%s=1.0_dp
4269	ENDIF
4270
4271	END SUBROUTINE resetpoly_R
4272
4273	SUBROUTINE resetpoly_RN(S2,FL,K)
4274	implicit none
4275	type (real_8),INTENT(INOUT),dimension(:)::S2
4276	logical(lp),INTENT(IN)::FL
4277
4278	INTEGER,optional,INTENT(IN)::k
4279	INTEGER J,i,N
4280
4281	if(present(k)) then
4282	I=LBOUND(S2,DIM=1)
4283	N=LBOUND(S2,DIM=1)+K-1
4284	else
4285	I=LBOUND(S2,DIM=1)
4286	N=UBOUND(S2,DIM=1)
4287	endif
4288
4289	DO J=I,N
4290	call resetpoly_R(s2(j),FL)
4291	enddo
4292
4293
4294	END SUBROUTINE resetpoly_RN
4295
4296	SUBROUTINE resetpoly_R31(S2) ! STAYS REAL FOR PTC
4297	implicit none
4298	integer ipause, mypauses
4299	type (real_8),INTENT(INOUT)::S2
4300
4301	if(s2%kind==3) then
4302	if(s2%alloc) then
4303	line= "Allocated in resetpoly_R31"
4304	ipause=mypauses(2,line)
4305	endif
4306	s2%kind=1
4307	s2%i=0
4308	s2%s=1.0_dp
4309	ENDIF
4310	END SUBROUTINE resetpoly_R31
4311
4312	SUBROUTINE resetpoly_R31N(S2,K)
4313	implicit none
4314	type (real_8),INTENT(INOUT),dimension(:)::S2
4315
4316	INTEGER,optional,INTENT(IN)::k
4317	INTEGER J,i,N
4318
4319	if(present(k)) then
4320	I=LBOUND(S2,DIM=1)
4321	N=LBOUND(S2,DIM=1)+K-1
4322	else
4323	I=LBOUND(S2,DIM=1)
4324	N=UBOUND(S2,DIM=1)
4325	endif
4326
4327	DO J=I,N
4328	call resetpoly_R31(s2(j))
4329	enddo
4330
4331
4332	END SUBROUTINE resetpoly_R31N
4333
4334
4335	SUBROUTINE resetpoly0(S2)
4336	implicit none
4337	type (real_8),INTENT(INOUT)::S2
4338
4339	if(s2%alloc) call killtpsa(s2%t)
4340	s2%alloc=f
4341	s2%kind=0
4342	s2%r=0.0_dp
4343
4344	s2%i=0
4345	s2%s=1.0_dp
4346
4347	END SUBROUTINE resetpoly0
4348
4349	SUBROUTINE K_OPT(S1,S2,s3,s4,s5,s6,s7,s8,s9,s10)
4350	implicit none
4351	type (real_8),INTENT(INout)::S1
4352	type (real_8),optional, INTENT(INout):: S2,s3,s4,s5,s6,s7,s8,s9,s10
4353	call resetpoly0(s1)
4354	if(present(s2)) call resetpoly0(s2)
4355	if(present(s3)) call resetpoly0(s3)
4356	if(present(s4)) call resetpoly0(s4)
4357	if(present(s5)) call resetpoly0(s5)
4358	if(present(s6)) call resetpoly0(s6)
4359	if(present(s7)) call resetpoly0(s7)
4360	if(present(s8)) call resetpoly0(s8)
4361	if(present(s9)) call resetpoly0(s9)
4362	if(present(s10))call resetpoly0(s10)
4363	END SUBROUTINE K_OPT
4364
4365
4366	SUBROUTINE Ke_OPT(S1,S2,s3,s4,s5,s6,s7,s8,s9,s10)
4367	implicit none
4368	type (env_8),INTENT(INout)::S1
4369	type (env_8),optional, INTENT(INout):: S2,s3,s4,s5,s6,s7,s8,s9,s10
4370	call resetenv(s1)
4371	if(present(s2)) call resetenv(s2)
4372	if(present(s3)) call resetenv(s3)
4373	if(present(s4)) call resetenv(s4)
4374	if(present(s5)) call resetenv(s5)
4375	if(present(s6)) call resetenv(s6)
4376	if(present(s7)) call resetenv(s7)
4377	if(present(s8)) call resetenv(s8)
4378	if(present(s9)) call resetenv(s9)
4379	if(present(s10))call resetenv(s10)
4380	END SUBROUTINE Ke_OPT
4381
4382
4383
4384
4385
4386
4387	! SUBROUTINE kill_c(k)
4388	! implicit none
4389	! integer ipause, mypauses
4390	! integer k
4391
4392	! call dallsta(exi) ! Etienne checking routine
4393	! if(exi/=ent) then
4394	! write(line,'(3(1x,i8))') ent,exi,k
4395	! ipause=mypauses(1999,line)
4396	! k=k*10000
4397	! endif
4398
4399	! END SUBROUTINE kill_c
4400
4401	SUBROUTINE resetpolyn0(S2,K)
4402	implicit none
4403	type (real_8),INTENT(INOUT),dimension(:)::S2
4404	INTEGER,optional,INTENT(IN)::k
4405	INTEGER J,i,N
4406
4407	if(present(k)) then
4408	I=LBOUND(S2,DIM=1)
4409	N=LBOUND(S2,DIM=1)+K-1
4410	else
4411	I=LBOUND(S2,DIM=1)
4412	N=UBOUND(S2,DIM=1)
4413	endif
4414
4415	DO J=I,N
4416	call resetpoly0(s2(j))
4417	enddo
4418
4419
4420	END SUBROUTINE resetpolyn0
4421
4422
4423	SUBROUTINE allocpoly(S2)
4424	implicit none
4425	type (real_8),INTENT(INOUT)::S2
4426
4427	!if(s2%alloc) call killtpsa(s2%t) ADDED ETIENNE
4428	s2%alloc=f
4429	s2%kind=1
4430	s2%r=0.0_dp
4431	s2%i=0
4432	s2%g=0
4433	s2%nb=0
4434	s2%s=1.0_dp
4435
4436	END SUBROUTINE allocpoly
4437
4438	SUBROUTINE A_OPT(S1,S2,s3,s4,s5,s6,s7,s8,s9,s10)
4439	implicit none
4440	type (real_8),INTENT(INout)::S1
4441	type (real_8),optional, INTENT(INout):: S2,s3,s4,s5,s6,s7,s8,s9,s10
4442	call allocpoly(s1)
4443	if(present(s2)) call allocpoly(s2)
4444	if(present(s3)) call allocpoly(s3)
4445	if(present(s4)) call allocpoly(s4)
4446	if(present(s5)) call allocpoly(s5)
4447	if(present(s6)) call allocpoly(s6)
4448	if(present(s7)) call allocpoly(s7)
4449	if(present(s8)) call allocpoly(s8)
4450	if(present(s9)) call allocpoly(s9)
4451	if(present(s10))call allocpoly(s10)
4452	END SUBROUTINE A_opt
4453
4454	! SUBROUTINE alloc_c
4455	! implicit none
4456	! Etienne checking routine
4457	! call dallsta(ent)
4458
4459
4460	! END SUBROUTINE alloc_c
4461
4462	SUBROUTINE allocpolyn(S2,K)
4463	implicit none
4464	type (real_8),INTENT(INOUT),dimension(:)::S2
4465	INTEGER,optional,INTENT(IN)::k
4466	INTEGER J,i,N
4467
4468	if(present(k)) then
4469	I=LBOUND(S2,DIM=1)
4470	N=LBOUND(S2,DIM=1)+K-1
4471	else
4472	I=LBOUND(S2,DIM=1)
4473	N=UBOUND(S2,DIM=1)
4474	endif
4475
4476	DO J=I,N
4477	call allocpoly(s2(j))
4478	enddo
4479
4480
4481	END SUBROUTINE allocpolyn
4482
4483
4484	subroutine init_map_p(NO1,ND1,NP1,NDPT1,log)
4485	implicit none
4486	integer NO1,ND1,NP1,NDPT1
4487	logical(lp) log
4488	call init_map_c(NO1,ND1,NP1,NDPT1,log)
4489	call set_in_polyp(log)
4490	end subroutine init_map_p
4491
4492	subroutine init_tpsa_p(NO1,NP1,log)
4493	implicit none
4494	integer NO1,NP1
4495	logical(lp) log
4496	call init_tpsa_c(NO1,NP1,log)
4497	call set_in_polyp(log)
4498	end subroutine init_tpsa_p
4499
4500
4501	subroutine set_in_polyp(log)
4502	implicit none
4503	logical(lp) log
4504	integer iia(4),icoast(4)
4505	call liepeek(iia,icoast)
4506	old=log
4507	NO=iia(1)
4508	ND=iia(3)
4509	ND2=iia(3)*2
4510	NP=iia(2)-nd2
4511	NDPT=icoast(4)
4512	NV=iia(2)
4513	! i_ =cmplx(zero,one,kind=dp)
4514	end subroutine set_in_polyp
4515
4516	SUBROUTINE EQUAL2D(S2,S1)
4517	implicit none
4518	type (real_8),INTENT(inOUT)::S2(:,:)
4519	type (real_8),INTENT(IN)::S1(:,:)
4520	integer i,J,I1(2),I2(2)
4521
4522	I1(1)=LBOUND(S1,DIM=1)
4523	I1(2)=LBOUND(S1,DIM=2)
4524	I2(1)=LBOUND(S2,DIM=1)
4525	I2(2)=LBOUND(S2,DIM=2)
4526	do i=I1(1),UBOUND(S1,DIM=1)
4527	do J=I1(2),UBOUND(S1,DIM=2)
4528	S2(I-I1(1)+I2(1),J-I1(2)+I2(2))=S1(I,J)
4529	ENDDO
4530	ENDDO
4531
4532	end SUBROUTINE EQUAL2D
4533
4534	SUBROUTINE EQUAL1D(S2,S1)
4535	implicit none
4536	type (real_8),INTENT(inOUT)::S2(:)
4537	type (real_8),INTENT(IN)::S1(:)
4538	integer i,I1,I2
4539
4540	I1=LBOUND(S1,DIM=1)
4541	I2=LBOUND(S2,DIM=1)
4542	do i=I1,UBOUND(S1,DIM=1)
4543	S2(I-I1+I2)=S1(I)
4544	ENDDO
4545
4546	end SUBROUTINE EQUAL1D
4547
4548	SUBROUTINE EQUAL(S2,S1)
4549	implicit none
4550	integer ipause, mypauses
4551	type (real_8),INTENT(inOUT)::S2
4552	type (real_8),INTENT(IN)::S1
4553	! integer localmaster
4554
4555
4556	if(s1%kind==0) then
4557	line= " You are putting kind=0 into something"
4558	ipause=mypauses(10,line)
4559	endif
4560	if(s2%kind==3.and.(.not.setknob)) then
4561	line= " You are putting something into a knob kind=3"
4562	ipause=mypauses(11,line)
4563	endif
4564
4565	if (s2%kind>0) then ! S2 exist
4566	if(S2%kind==S1%kind) then
4567	select case(S1%kind)
4568	case(m1)
4569	! localmaster=master ! added because of knob problems 2001.9.19
4570	! call check_snake ! added because of knob problems 2001.9.19
4571	! master=0
4572	S2%R=S1%R
4573	! master=localmaster ! added because of knob problems 2001.9.19
4574	case(m2)
4575	! localmaster=master
4576	call check_snake
4577	! 2002.12.26 master=0
4578	S2%t=S1%t
4579	! master=localmaster
4580	case(m3)
4581	s2%r=S1%r ! Knob stays a knob and real stays real 2002.10.9
4582	! line= " You are putting kind=3 (TPSA) into another kind=3"
4583	! ipause=mypauses(13,line)
4584	end select
4585	elseif(S2%kind>S1%kind ) then
4586	if(S1%kind/=2) then
4587	s2%r=S1%r
4588	if(s2%kind/=3) s2%kind=1 !Knob stays a knob and real stays real 2002.10.9
4589	else
4590	s2%r=S1%t.sub.'0' ! setting a kind=3
4591	endif
4592	else
4593	select case(s1%kind)
4594	case(m2)
4595	if(.not.s2%alloc) then
4596	call alloc(s2%t)
4597	s2%alloc=t
4598	endif
4599	s2%kind=2
4600	! localmaster=master
4601	call check_snake
4602	! 2002.12.26 master=0
4603
4604	S2%t=S1%t
4605
4606	! master=localmaster
4607	case(m3)
4608	if(.not.s2%alloc) then
4609	call alloc(s2%t)
4610	s2%alloc=t
4611	endif
4612	s2%kind=2
4613	! localmaster=master
4614	call check_snake
4615	! 2002.12.26 master=0
4616	if(knob) then
4617	call varfk1(S1)
4618	S2%t=varf1
4619	else
4620	S2%R=S1%R
4621	s2%kind=1
4622	endif
4623	! master=localmaster
4624	end select
4625	endif
4626
4627	else ! S2 does not exist
4628
4629
4630
4631	if(S1%kind==1) then ! what is s1
4632	if(s2%i==0) then
4633	S2%R=S1%R
4634	s2%kind=1
4635	elseif(s2%i>0.and.s2%i<=nv) then
4636	call alloc(s2%t)
4637	s2%t=(/S1%R,S2%S/).var.s2%i
4638	! call var(s2%t,S1%R,S2%S,s2%i)
4639	! s2%i=0
4640	s2%kind=2
4641	s2%alloc=t
4642	else
4643	write(6,*) "EQUAL IN m_POLYMORPH ",s2%i,S2%KIND,S2%ALLOC
4644	WRITE(6,*) " I "
4645	READ(5,*) s2%i
4646	WRITE(6,*) SQRT(FLOAT(s2%i))
4647	stop 777
4648	endif
4649
4650	else
4651	if(insane_PTC) then
4652	if(.not.s2%alloc) then
4653	call alloc(s2%t)
4654	endif
4655	! localmaster=master
4656	call check_snake
4657	! 2002.12.26 master=0
4658	S2%t=S1%t
4659	! master=localmaster
4660	s2%kind=2
4661	s2%alloc=t
4662	else
4663	w_p=0
4664	w_p%nc=5
4665	w_p%fc='(4(1X,A72,/),(1X,A72))'
4666	write(w_p%c(1),'(A23,I4,A19)') " You are putting kind= ", s1%kind," (TPSA) in a kind=0"
4667	w_p%c(2)= " We do not allow that anymore for safety reasons"
4668	w_p%c(3)= " If you insist on it, modify real_polymorph and complex_polymorph"
4669	w_p%c(4)= " at your own insane risk "
4670	w_p%c(5)= " Etienne Forest/Frank Schmidt"
4671	! call !write_e(778)
4672	w_p%nc=sqrt(-float(w_p%nc))
4673	endif
4674	endif ! end of what is s1
4675
4676	endif ! S2 does not exist
4677
4678	END SUBROUTINE EQUAL
4679
4680	SUBROUTINE complexreal_8(S2,S1)
4681	implicit none
4682	complex(dp) ,INTENT(inOUT)::S2
4683	type (real_8),INTENT(IN)::S1
4684	! integer localmaster
4685
4686
4687	select case(S1%kind)
4688	case(m1)
4689	S2=S1%R
4690	case(m2)
4691	! localmaster=master
4692	call check_snake
4693	! 2002.12.26 master=0
4694	S2=S1%t.sub.'0'
4695	! master=localmaster
4696	case(m3)
4697	! localmaster=master
4698	call check_snake
4699	! 2002.12.26 master=0
4700	S2=S1%r
4701	! master=localmaster
4702	case default
4703	w_p=0
4704	w_p%nc=2
4705	w_p%fc='((1X,A72,/,1x,a72))'
4706	w_p%fi='(2((1X,i4)))'
4707	w_p%c(1)= " trouble in complexreal_8 "
4708	w_p%c(2)= "s1%kind "
4709	w_p=(/s1%kind /)
4710	! call !write_e(0)
4711	end select
4712
4713	END SUBROUTINE complexreal_8
4714
4715
4716
4717	SUBROUTINE realEQUAL(S2,S1)
4718	implicit none
4719	real(dp),INTENT(inOUT)::S2
4720	type (real_8),INTENT(IN)::S1
4721	! integer localmaster
4722
4723
4724	select case(S1%kind)
4725	case(m1)
4726	S2=S1%R
4727	case(m2)
4728	! localmaster=master
4729	call check_snake
4730	! 2002.12.26 master=0
4731	S2=S1%t.sub.'0'
4732	! master=localmaster
4733	case(m3)
4734	! localmaster=master
4735	call check_snake
4736	! 2002.12.26 master=0
4737	S2=S1%r
4738	! master=localmaster
4739	case default
4740	w_p=0
4741	w_p%nc=2
4742	w_p%fc='((1X,A72,/,1x,a72))'
4743	w_p%fi='(2((1X,i4)))'
4744	w_p%c(1)= " trouble in realEQUAL "
4745	w_p%c(2)= "s1%kind "
4746	w_p=(/s1%kind /)
4747	! call !write_e(0)
4748	end select
4749
4750	END SUBROUTINE realEQUAL
4751
4752	SUBROUTINE singleEQUAL(S2,S1)
4753	implicit none
4754	real(sp),INTENT(inOUT)::S2
4755	type (real_8),INTENT(IN)::S1
4756	! integer localmaster
4757
4758
4759	select case(S1%kind)
4760	case(m1)
4761	S2=S1%R
4762	case(m2)
4763	! localmaster=master
4764	call check_snake
4765	! 2002.12.26 master=0
4766	S2=S1%t.sub.'0'
4767	! master=localmaster
4768	case(m3)
4769	! localmaster=master
4770	call check_snake
4771	! 2002.12.26 master=0
4772	S2=S1%r
4773	! master=localmaster
4774	case default
4775	w_p=0
4776	w_p%nc=2
4777	w_p%fc='((1X,A72,/,1x,a72))'
4778	w_p%fi='(2((1X,i4)))'
4779	w_p%c(1)= " trouble in realEQUAL "
4780	w_p%c(2)= "s1%kind "
4781	w_p=(/s1%kind /)
4782	! call !write_e(0)
4783	end select
4784
4785	END SUBROUTINE singleEQUAL
4786
4787
4788
4789	SUBROUTINE taylorEQUAL(S2,S1)
4790	implicit none
4791	type (taylor),INTENT(inOUT)::S2
4792	type (real_8),INTENT(IN)::S1
4793	! integer localmaster
4794
4795
4796	select case(S1%kind)
4797	case(m1)
4798	S2=S1%R
4799	case(m2)
4800	! localmaster=master
4801	call check_snake
4802	! 2002.12.26 master=0
4803	S2=S1%t
4804	! master=localmaster
4805	case(m3)
4806	! localmaster=master
4807	call check_snake
4808	! 2002.12.26 master=0
4809	if(knob) then
4810	call varfk1(S1)
4811	S2=varf1
4812	else
4813	S2=S1%R
4814	endif
4815	! master=localmaster
4816	case default
4817	w_p=0
4818	w_p%nc=2
4819	w_p%fc='((1X,A72,/,1x,a72))'
4820	w_p%fi='(2((1X,i4)))'
4821	w_p%c(1)= " trouble in taylorEQUAL "
4822	w_p%c(2)= "s1%kind "
4823	w_p=(/s1%kind /)
4824	! call !write_e(0)
4825	end select
4826
4827	END SUBROUTINE taylorEQUAL
4828
4829	SUBROUTINE EQUALtaylor(S2,S1)
4830	implicit none
4831	integer ipause, mypauses
4832	type (real_8),INTENT(inOUT)::S2
4833	type (taylor),INTENT(IN)::S1
4834	! integer localmaster
4835
4836	IF(S2%KIND==3.and.(.not.setknob)) THEN ! 2002.10.9
4837	line= "Forbidden in EQUALtaylor: s2 is a knob "
4838	ipause=mypauses(0,line)
4839	ENDIF
4840
4841	! localmaster=master
4842	call check_snake
4843	! 2002.12.26 master=0
4844	if(s2%kind/=3) then
4845	if(.not.s2%alloc) then
4846	call alloc(s2%t)
4847	s2%alloc=t
4848	endif
4849	S2%t=S1
4850	s2%KIND=2
4851	else
4852	s2%r=S1.sub.'0' ! 2002.10.9
4853	endif
4854	! master=localmaster
4855
4856	END SUBROUTINE EQUALtaylor
4857
4858	SUBROUTINE real_8univ(S2,S1)
4859	implicit none
4860	integer ipause, mypauses
4861	type (real_8),INTENT(inOUT)::S2
4862	type (universal_taylor),INTENT(IN)::S1
4863	! integer localmaster
4864
4865	IF(S2%KIND==M3) THEN
4866	line= "Forbidden in real_8univ: s2 is a knob "
4867	ipause=mypauses(0,line)
4868	ENDIF
4869
4870	! localmaster=master
4871	! call check_snake
4872	! 2002.12.26 master=0
4873	if(.not.s2%alloc) then
4874	call alloc(s2%t)
4875	s2%alloc=t
4876	endif
4877	S2%t=S1
4878	s2%KIND=2
4879
4880	! master=localmaster
4881	END SUBROUTINE real_8univ
4882
4883	SUBROUTINE univreal_8(S2,S1) ! new sagan
4884	implicit none
4885	type (universal_taylor),INTENT(inOUT)::S2
4886	type (real_8),INTENT(IN)::S1
4887	! integer localmaster
4888
4889	select case(S1%kind)
4890	case(m1)
4891	S2=S1%R
4892	case(m2)
4893	! localmaster=master
4894	call check_snake
4895	! 2002.12.26 master=0
4896	S2=S1%t
4897	! master=localmaster
4898	case(m3)
4899	! localmaster=master
4900	call check_snake
4901	! 2002.12.26 master=0
4902	if(knob) then
4903	call varfk1(S1)
4904	S2=varf1
4905	else
4906	S2=S1%R
4907	endif
4908	! master=localmaster
4909	case default
4910	w_p=0
4911	w_p%nc=2
4912	w_p%fc='((1X,A72,/,1x,a72))'
4913	w_p%fi='(2((1X,i4)))'
4914	w_p%c(1)= " trouble in univreal_8 "
4915	w_p%c(2)= "s1%kind "
4916	w_p=(/s1%kind /)
4917	! call !write_e(0)
4918	end select
4919	END SUBROUTINE univreal_8
4920
4921
4922
4923	SUBROUTINE Dequaldacon(S2,R1)
4924	implicit none
4925	integer ipause, mypauses
4926	type (real_8),INTENT(inOUT)::S2
4927	real(dp),INTENT(IN)::R1
4928
4929	IF(S2%KIND==M3) THEN
4930	if(setknob) then
4931	s2%r=r1
4932	return
4933	else
4934	line= "Forbidden in Dequaldacon: s2 is a knob "
4935	ipause=mypauses(0,line)
4936	endif
4937	ENDIF
4938
4939
4940	!localmaster=master
4941	! call check_snake
4942	! master=0
4943	if (s2%kind/=0) then ! S2 exist
4944	if(S2%kind==1.or.s2%kind==3) then
4945	S2%R=R1
4946	else
4947	! S2%t=R1
4948	s2%r=r1
4949	s2%kind=1
4950	! call kill(s2%t)
4951	! s2%alloc=f
4952	endif
4953
4954	else ! S2 does not exist
4955
4956
4957
4958	if(s2%i==0) then
4959	S2%R=R1
4960	s2%kind=1
4961	elseif(s2%i>0.and.s2%i<=nv) then
4962	call alloc(s2%t)
4963	s2%t=(/R1,s2%S/).var.s2%i
4964	! call var(s2%t,R1,s2%S,s2%i)
4965	! s2%i=0
4966	s2%kind=2
4967	s2%alloc=t
4968	else
4969	line= "trouble in Dequaldacon in Real_polymorph"
4970	ipause=mypauses(-788,line)
4971	endif
4972
4973
4974	endif ! S2 not allocated
4975	! master=localmaster
4976	END SUBROUTINE Dequaldacon
4977
4978	SUBROUTINE iequaldaconn(S2,R1)
4979	implicit none
4980	type (real_8),INTENT(inOUT),dimension(:)::S2
4981	integer,INTENT(IN)::R1
4982	integer i
4983	! integer localmaster
4984
4985	! localmaster=master
4986	! call check_snake
4987	! 2002.12.26 master=0
4988	do i=1,r1
4989	call resetpoly(s2(i))
4990	s2(i)%i=i
4991	enddo
4992	! master=localmaster
4993	END SUBROUTINE iequaldaconn
4994
4995
4996	SUBROUTINE equaldacon(S2,R1)
4997	implicit none
4998	integer ipause, mypauses
4999	type (real_8),INTENT(inOUT)::S2
5000	real(sp),INTENT(IN)::R1
5001
5002	if(real_warning) call real_stop
5003	IF(S2%KIND==M3) THEN
5004	if(setknob) then
5005	s2%r=REAL(R1,kind=DP)
5006	return
5007	else
5008	line= "Forbidden in equaldacon: s2 is a knob "
5009	ipause=mypauses(0,line)
5010	endif
5011	ENDIF
5012	!localmaster=master
5013	! call check_snake
5014	! master=0
5015	if (s2%kind/=0) then ! S2 exist
5016	if(S2%kind==1.or.s2%kind==3) then
5017	S2%R=REAL(R1,kind=DP)
5018	else
5019	! S2%t=REAL(R1,kind=DP)
5020	! S2%t=R1
5021	s2%r=REAL(R1,kind=DP)
5022	s2%kind=1
5023	! call kill(s2%t)
5024	! s2%alloc=f
5025	endif
5026
5027	else ! S2 does not exist
5028
5029
5030
5031	if(s2%i==0) then
5032	S2%R=REAL(R1,kind=DP)
5033	s2%kind=1
5034	elseif(s2%i>0.and.s2%i<=nv) then
5035	call alloc(s2%t)
5036	s2%t=(/REAL(R1,kind=DP),S2%S/).var.s2%i
5037	! call var(s2%t,REAL(R1,kind=DP),S2%S,s2%i)
5038	! s2%i=0
5039	s2%kind=2
5040	s2%alloc=t
5041	else
5042	stop 779
5043	endif
5044
5045
5046	endif ! S2 not allocated
5047	! master=localmaster
5048	END SUBROUTINE equaldacon
5049
5050	SUBROUTINE iequaldacon(S2,R1)
5051	implicit none
5052	integer ipause, mypauses
5053	type (real_8),INTENT(inOUT)::S2
5054	integer,INTENT(IN)::R1
5055
5056	IF(S2%KIND==M3) THEN
5057	if(setknob) then
5058	s2%r=r1
5059	return
5060	else
5061	line= "Forbidden in iequaldacon: s2 is a knob "
5062	ipause=mypauses(0,line)
5063	endif
5064	ENDIF
5065	!localmaster=master
5066	! call check_snake
5067	! master=0
5068	if (s2%kind/=0) then ! S2 exist
5069	if(S2%kind==1.or.s2%kind==3) then
5070	S2%R=REAL(R1,kind=DP)
5071	else
5072	! S2%t=REAL(R1,kind=DP)
5073	s2%r=REAL(R1,kind=DP)
5074	s2%kind=1
5075	! call kill(s2%t)
5076	! s2%alloc=f
5077	endif
5078
5079	else ! S2 does not exist
5080
5081
5082
5083	if(s2%i==0) then
5084	S2%R=REAL(R1,kind=DP)
5085	s2%kind=1
5086	elseif(s2%i>0.and.s2%i<=nv) then
5087	call alloc(s2%t)
5088	s2%t=(/REAL(R1,kind=DP),S2%S/).var.s2%i
5089	! call var(s2%t,REAL(R1,kind=DP),S2%S,s2%i)
5090	! s2%i=0
5091	s2%kind=2
5092	s2%alloc=t
5093	else
5094	stop 776
5095	endif
5096
5097
5098	endif ! S2 not allocated
5099	! master=localmaster
5100	END SUBROUTINE iequaldacon
5101
5102	subroutine assp(s1)
5103	implicit none
5104	TYPE (real_8) s1
5105	integer ipause,mypauses
5106	! lastmaster=master ! 2002.12.13
5107
5108	select case(master)
5109	case(0:ndumt-1)
5110	master=master+1
5111	case(ndumt)
5112	line= " cannot indent anymore "
5113	ipause=mypauses(0,line)
5114	end select
5115	! write(26,*) " real polymorph ",master
5116
5117	call ass0(s1%t)
5118	s1%alloc=t
5119	s1%kind=2
5120	s1%i=0
5121
5122	end subroutine ASSp
5123
5124	subroutine assp_no_master(s1)
5125	implicit none
5126	TYPE (real_8) s1
5127
5128
5129	call ass0(s1%t)
5130	s1%alloc=t
5131	s1%kind=1
5132	s1%i=0
5133
5134	end subroutine assp_no_master
5135
5136	subroutine assp_master(s1)
5137	implicit none
5138	TYPE (real_8) s1
5139	integer ipause,mypauses
5140	! lastmaster=master ! 2002.12.13
5141
5142	select case(master)
5143	case(0:ndumt-1)
5144	master=master+1
5145	case(ndumt)
5146	line= " cannot indent anymore "
5147	ipause=mypauses(0,line)
5148	end select
5149	! write(26,*) " real polymorph ",master
5150
5151	call ass0(s1%t)
5152	s1%alloc=t
5153	s1%kind=1
5154	s1%i=0
5155
5156	end subroutine assp_master
5157
5158	! complex
5159
5160	FUNCTION datant( S1 )
5161	implicit none
5162	TYPE (real_8) datant
5163	TYPE (real_8), INTENT (IN) :: S1
5164	TYPE (complextaylor) w
5165	integer localmaster
5166
5167	select case(s1%kind)
5168	case(m1)
5169	datant%r=atan(s1%r)
5170	datant%kind=1
5171	case(m2)
5172	localmaster=master
5173	call ass(datant)
5174	call alloc(w)
5175	w%r=s1%t
5176	w= atan(w)
5177	datant%t=w%r
5178	call kill(w)
5179	master=localmaster
5180	case(m3)
5181	if(knob) then
5182	localmaster=master
5183	call ass(datant)
5184	call alloc(w)
5185	call varfk1(S1)
5186	w%r=varf1
5187	w= atan(w)
5188	datant%t=w%r
5189	call kill(w)
5190	master=localmaster
5191	else
5192	datant%r=atan(s1%r)
5193	datant%kind=1
5194	endif
5195	case default
5196	w_p=0
5197	w_p%nc=2
5198	w_p%fc='((1X,A72,/,1x,a72))'
5199	w_p%fi='(2((1X,i4)))'
5200	w_p%c(1)= " trouble in datant "
5201	w_p%c(2)= "s1%kind "
5202	w_p=(/s1%kind /)
5203	! call !write_e(0)
5204	end select
5205	END FUNCTION datant
5206
5207	FUNCTION datanht( S1 )
5208	implicit none
5209	TYPE (real_8) datanht
5210	TYPE (real_8), INTENT (IN) :: S1
5211	integer localmaster
5212
5213	select case(s1%kind)
5214	case(m1)
5215	datanht%r=log((1+s1%r)/(1-s1%r))/2.0_dp
5216	datanht%kind=1
5217	case(m2)
5218	localmaster=master
5219	call ass(datanht)
5220	datanht%t=atanh(s1%t)
5221	master=localmaster
5222	case(m3)
5223	if(knob) then
5224	localmaster=master
5225	call ass(datanht)
5226	call varfk1(S1)
5227	datanht%t=atanh(varf1)
5228	master=localmaster
5229	else
5230	datanht%r=log((1+s1%r)/sqrt(1-s1%r))/2.0_dp
5231	datanht%kind=1
5232	endif
5233	case default
5234	w_p=0
5235	w_p%nc=2
5236	w_p%fc='((1X,A72,/,1x,a72))'
5237	w_p%fi='(2((1X,i4)))'
5238	w_p%c(1)= " trouble in datant "
5239	w_p%c(2)= "s1%kind "
5240	w_p=(/s1%kind /)
5241	! call !write_e(0)
5242	end select
5243	END FUNCTION datanht
5244
5245
5246
5247
5248
5249	FUNCTION datanDt( S1 )
5250	implicit none
5251	TYPE (real_8) datanDt
5252	TYPE (real_8), INTENT (IN) :: S1
5253	TYPE (complextaylor) w
5254	integer localmaster
5255
5256	select case(s1%kind)
5257	case(m1)
5258	datanDt%r=atan(s1%r)*RAD_TO_DEG_
5259	datanDt%kind=1
5260	case(m2)
5261	localmaster=master
5262	call ass(datanDt)
5263	call alloc(w)
5264	w%r=s1%t
5265	w= atan(w)
5266	datanDt%t=w%r*RAD_TO_DEG_
5267	call kill(w)
5268	master=localmaster
5269	case(m3)
5270	if(knob) then
5271	localmaster=master
5272	call ass(datanDt)
5273	call alloc(w)
5274	call varfk1(S1)
5275	w%r=varf1
5276	w= atan(w)
5277	datanDt%t=w%r*RAD_TO_DEG_
5278	call kill(w)
5279	master=localmaster
5280	else
5281	datanDt%r=atan(s1%r)*RAD_TO_DEG_
5282	datanDt%kind=1
5283	endif
5284	case default
5285	w_p=0
5286	w_p%nc=2
5287	w_p%fc='((1X,A72,/,1x,a72))'
5288	w_p%fi='(2((1X,i4)))'
5289	w_p%c(1)= " trouble in datanDt "
5290	w_p%c(2)= "s1%kind "
5291	w_p=(/s1%kind /)
5292	! call !write_e(0)
5293	end select
5294	END FUNCTION datanDt
5295
5296	FUNCTION absoftdatanDr( S1 )
5297	implicit none
5298	real(dp) absoftdatanDr
5299	real(dp), INTENT (IN) :: S1
5300	absoftdatanDr=atan(s1)*RAD_TO_DEG_
5301	END FUNCTION absoftdatanDr
5302
5303	FUNCTION datan2t( S2,S1 )
5304	implicit none
5305	integer ipause, mypauses
5306	TYPE (real_8) datan2t
5307	TYPE (real_8), INTENT (IN) :: S2,S1
5308	TYPE (complextaylor) w
5309	real(dp) ANG
5310	integer localmaster ,si
5311	si=1.0_dp
5312	if(s2%kind==0.or.S1%kind==0) then
5313	line= " Problems in datan2t "
5314	ipause=mypauses(0,line)
5315	endif
5316	if(S2%kind==S1%kind) then ! 1
5317	select case(S1%kind)
5318	case(m1)
5319	DATAN2T%R=ATAN2(S2%R,S1%R)
5320	datan2t%kind=1
5321	case(m2)
5322	localmaster=master
5323	call ass(datan2t)
5324	call alloc(w)
5325	ANG=ATAN2(S2%T.SUB.'0',S1%T.SUB.'0')
5326	if(pil>abs(ANG).or.abs(ANG)>pim) then
5327	w%r=s2%t/s1%t
5328	w=atan(w)
5329	datan2t%t=w%r
5330	datan2t%t=datan2t%t-(datan2t%t.SUB.'0')+ ANG
5331	else
5332	if((s2%t.sub.'0')<0.0_dp) si=-1.0_dp
5333	w%r=s1%t/SQRT(s2%t2+s1%t2)
5334	w=acos(w)
5335	datan2t%t=w%r
5336	datan2t%t=si*(datan2t%t-(datan2t%t.SUB.'0'))+ ANG
5337	endif
5338	call kill(w)
5339	master=localmaster
5340	case(m3)
5341	localmaster=master
5342	if(knob) then !knob 1
5343	call ass(datan2t)
5344	call alloc(w)
5345	call varfk1(S1)
5346	call varfk2(S2)
5347	ANG=ATAN2(varf2.SUB.'0',varf1.SUB.'0')
5348	if(pil>abs(ANG).or.abs(ANG)>pim) then
5349	w%r=varf2/varf1
5350	w=atan(w)
5351	datan2t%t=w%r
5352	datan2t%t=datan2t%t-(datan2t%t.SUB.'0')+ ANG
5353	else
5354	if((varf2.sub.'0')<0.0_dp) si=-1.0_dp
5355	w%r=varf1/SQRT(varf22+varf12)
5356	w=acos(w)
5357	datan2t%t=w%r
5358	datan2t%t=si*(datan2t%t-(datan2t%t.SUB.'0'))+ ANG
5359	endif
5360	call kill(w)
5361	else !knob 1
5362	DATAN2T%R=ATAN2(S2%R,S1%R)
5363	datan2t%kind=1
5364	endif !knob 1
5365	master=localmaster
5366	end select
5367	elseif(S2%kind>S1%kind) then !1
5368	localmaster=master
5369	! call ass(datan2t)
5370	call alloc(w)
5371	if(s2%kind==2) then !
5372	call ass(datan2t)
5373	ANG=ATAN2(S2%T.SUB.'0',S1%R)
5374	if(pil>abs(ANG).or.abs(ANG)>pim) then
5375	w%r=s2%t/s1%r
5376	w=atan(w)
5377	datan2t%t=w%r
5378	datan2t%t=datan2t%t-(datan2t%t.SUB.'0')+ ANG
5379
5380	else
5381	if((s2%t.sub.'0')<0.0_dp) si=-1.0_dp
5382	w%r=s1%t/SQRT(s2%t2+s1%r2)
5383	w=acos(w)
5384	datan2t%t=w%r
5385	datan2t%t=si*(datan2t%t-(datan2t%t.SUB.'0'))+ ANG
5386	endif
5387	else !
5388	if(knob) then ! knob
5389	call ass(datan2t)
5390	call varfk2(S2)
5391	if(s1%kind==1) then !!
5392	ANG=ATAN2(varf2.SUB.'0',S1%R)
5393	if(pil>abs(ANG).or.abs(ANG)>pim) then !!!
5394	w%r=varf2/s1%r
5395	w=atan(w)
5396	datan2t%t=w%r
5397	datan2t%t=datan2t%t-(datan2t%t.SUB.'0')+ ANG
5398
5399	else !!!
5400	if((varf2.sub.'0')<0.0_dp) si=-1.0_dp
5401	w%r=s1%R/SQRT(varf22+s1%R2) ! etienne error here????
5402	w=acos(w)
5403	datan2t%t=w%r
5404	datan2t%t=si*(datan2t%t-(datan2t%t.SUB.'0'))+ ANG
5405	endif !!!
5406
5407	else !!
5408	ANG=ATAN2(varf2.SUB.'0',S1%t.sub.'0')
5409	if(pil>abs(ANG).or.abs(ANG)>pim) then !!!
5410	w%r=varf2/s1%t
5411	w=atan(w)
5412	datan2t%t=w%r
5413	datan2t%t=datan2t%t-(datan2t%t.SUB.'0')+ ANG
5414
5415	else !!!
5416	if((varf2.sub.'0')<0.0_dp) si=-1.0_dp
5417	w%r=s1%t/SQRT(varf22+s1%t2)
5418	w=acos(w)
5419	datan2t%t=w%r
5420	datan2t%t=si*(datan2t%t-(datan2t%t.SUB.'0'))+ ANG
5421	endif !!!
5422	endif !!
5423	else ! knob
5424
5425	if(s1%kind==1) then !!
5426	DATAN2T%R=ATAN2(S2%R,S1%R)
5427	datan2t%kind=1
5428
5429	else !!
5430	call ass(datan2t)
5431	ANG=ATAN2(S2%R,S1%t.sub.'0')
5432	if(pil>abs(ANG).or.abs(ANG)>pim) then !!!
5433	w%r=S2%R/s1%t
5434	w=atan(w)
5435	datan2t%t=w%r
5436	datan2t%t=datan2t%t-(datan2t%t.SUB.'0')+ ANG
5437
5438	else !!!
5439	if((S2%R)<0.0_dp) si=-1.0_dp
5440	w%r=s1%t/SQRT(S2%R2+s1%t2)
5441	w=acos(w)
5442	datan2t%t=w%r
5443	datan2t%t=si*(datan2t%t-(datan2t%t.SUB.'0'))+ ANG
5444	endif !!!
5445	endif !!
5446
5447
5448	endif ! knob
5449	endif !
5450
5451	call kill(w)
5452	master=localmaster
5453	ELSE ! 1
5454	localmaster=master
5455	! call ass(datan2t)
5456	call alloc(w)
5457	if(s1%kind==2) then !
5458	call ass(datan2t)
5459	ANG=ATAN2(S2%R,S1%T.SUB.'0')
5460	if(pil>abs(ANG).or.abs(ANG)>pim) then
5461	w%r=s2%r/s1%t
5462	w=atan(w)
5463	datan2t%t=w%r
5464	datan2t%t=datan2t%t-(datan2t%t.SUB.'0')+ ANG
5465	else
5466	if(s2%r<0.0_dp) si=-1.0_dp
5467	w%r=s1%t/SQRT(s2%r2+s1%t2)
5468	w=acos(w)
5469	datan2t%t=w%r
5470	datan2t%t=si*(datan2t%t-(datan2t%t.SUB.'0'))+ ANG
5471	endif
5472	else !
5473	IF(KNOB) THEN ! KNOB 2
5474	call ass(datan2t)
5475	if(s2%kind==1) then !!
5476	call varfk1(S1)
5477	ANG=ATAN2(S2%R,varf1.SUB.'0')
5478	if(pil>abs(ANG).or.abs(ANG)>pim) then !!!
5479	w%r=s2%r/varf1
5480	w=atan(w)
5481	datan2t%t=w%r
5482	datan2t%t=datan2t%t-(datan2t%t.SUB.'0')+ ANG
5483	else !!!
5484	if(s2%r<0.0_dp) si=-1.0_dp
5485	w%r=varf1/SQRT(s2%r2+varf12)
5486	w=acos(w)
5487	datan2t%t=w%r
5488	datan2t%t=si*(datan2t%t-(datan2t%t.SUB.'0'))+ ANG
5489	endif !!!
5490	else !! s2%kind ==2
5491	call varfk1(S1)
5492	ANG=ATAN2(S2%t.sub.'0',varf1.SUB.'0')
5493	if(pil>abs(ANG).or.abs(ANG)>pim) then !!!
5494	w%r=s2%t/varf1
5495	w=atan(w)
5496	datan2t%t=w%r
5497	datan2t%t=datan2t%t-(datan2t%t.SUB.'0')+ ANG
5498	else !!!
5499	if((S2%t.sub.'0')<0.0_dp) si=-1.0_dp
5500	w%r=varf1/SQRT(s2%t2+varf12)
5501	w=acos(w)
5502	datan2t%t=w%r
5503	datan2t%t=si*(datan2t%t-(datan2t%t.SUB.'0'))+ ANG
5504	endif !!!
5505	endif !!
5506
5507	ELSE !KNOB 2
5508	if(s2%kind==1) then !!
5509	DATAN2T%R=ATAN2(S2%R,S1%R)
5510	datan2t%kind=1
5511	else !! s2%kind ==2
5512	call ass(datan2t)
5513	ANG=ATAN2(S2%t.sub.'0',S1%r)
5514	if(pil>abs(ANG).or.abs(ANG)>pim) then !!!
5515	w%r=s2%t/S1%r
5516	w=atan(w)
5517	datan2t%t=w%r
5518	datan2t%t=datan2t%t-(datan2t%t.SUB.'0')+ ANG
5519	else !!!
5520	if((S2%t.sub.'0')<0.0_dp) si=-1.0_dp
5521	w%r=S1%r/SQRT(s2%t2+S1%r2)
5522	w=acos(w)
5523	datan2t%t=w%r
5524	datan2t%t=si*(datan2t%t-(datan2t%t.SUB.'0'))+ ANG
5525	endif !!!
5526	endif !!
5527
5528
5529
5530
5531	ENDIF ! KNOB 2
5532
5533
5534	endif !
5535
5536	call kill(w)
5537	master=localmaster
5538	endif
5539
5540
5541	END FUNCTION datan2t
5542
5543	FUNCTION datan2dt( S2,S1 )
5544	implicit none
5545	integer ipause, mypauses
5546	TYPE (real_8) datan2dt
5547	TYPE (real_8), INTENT (IN) :: S2,S1
5548	TYPE (complextaylor) w
5549	real(dp) ANG
5550	integer localmaster ,si
5551	si=1.0_dp
5552	if(s2%kind==0.or.S1%kind==0) then
5553	line= " Problems in datand2t "
5554	ipause=mypauses(0,line)
5555	endif
5556	if(S2%kind==S1%kind) then
5557	select case(S1%kind)
5558	case(m1)
5559	datan2dt%R=ATAN2(S2%R,S1%R)*RAD_TO_DEG_
5560	datan2dt%kind=1
5561	case(m2)
5562	localmaster=master
5563	call ass(datan2dt)
5564	call alloc(w)
5565	ANG=ATAN2(S2%T.SUB.'0',S1%T.SUB.'0')
5566	if(pil>abs(ANG).or.abs(ANG)>pim) then
5567	w%r=s2%t/s1%t
5568	w=atan(w)
5569	datan2dt%t=w%r
5570	datan2dt%t=datan2dt%t-(datan2dt%t.SUB.'0')+ ANG
5571	datan2dt%t=datan2dt%t*RAD_TO_DEG_
5572	else
5573	if((s2%t.sub.'0')<0.0_dp) si=-1.0_dp
5574	w%r=s1%t/SQRT(s2%t2+s1%t2)
5575	w=acos(w)
5576	datan2dt%t=w%r
5577	datan2dt%t=si*(datan2dt%t-(datan2dt%t.SUB.'0'))+ ANG
5578	datan2dt%t=datan2dt%t*RAD_TO_DEG_
5579	endif
5580	call kill(w)
5581	master=localmaster
5582	case(m3)
5583	localmaster=master
5584	if(knob) then ! knob 1
5585	call ass(datan2dt)
5586	call alloc(w)
5587	call varfk1(S1)
5588	call varfk2(S2)
5589	ANG=ATAN2(varf2.SUB.'0',varf1.SUB.'0') !etienne
5590	if(pil>abs(ANG).or.abs(ANG)>pim) then
5591	w%r=varf2/varf1
5592	w=atan(w)
5593	datan2dt%t=w%r
5594	datan2dt%t=datan2dt%t-(datan2dt%t.SUB.'0')+ ANG
5595	datan2dt%t=datan2dt%t*RAD_TO_DEG_
5596	else
5597	if((varf2.sub.'0')<0.0_dp) si=-1.0_dp
5598	w%r=varf1/SQRT(varf22+varf12)
5599	w=acos(w)
5600	datan2dt%t=w%r
5601	datan2dt%t=si*(datan2dt%t-(datan2dt%t.SUB.'0'))+ ANG
5602	datan2dt%t=datan2dt%t*RAD_TO_DEG_
5603	endif
5604	call kill(w)
5605
5606	else ! knob 1
5607
5608	datan2dt%R=ATAN2(S2%R,S1%R)*RAD_TO_DEG_
5609	datan2dt%kind=1
5610
5611	endif ! knob 1
5612
5613	master=localmaster
5614	end select
5615	elseif(S2%kind>S1%kind) then
5616	localmaster=master
5617	! call ass(datan2dt)
5618	call alloc(w)
5619	if(s2%kind==2) then !
5620	call ass(datan2dt)
5621	ANG=ATAN2(S2%T.SUB.'0',S1%R)
5622	if(pil>abs(ANG).or.abs(ANG)>pim) then
5623	w%r=s2%t/s1%r
5624	w=atan(w)
5625	datan2dt%t=w%r
5626	datan2dt%t=datan2dt%t-(datan2dt%t.SUB.'0')+ ANG
5627	datan2dt%t=datan2dt%t*RAD_TO_DEG_
5628	else
5629	if((s2%t.sub.'0')<0.0_dp) si=-1.0_dp
5630	w%r=s1%t/SQRT(s2%t2+s1%r2)
5631	w=acos(w)
5632	datan2dt%t=w%r
5633	datan2dt%t=si*(datan2dt%t-(datan2dt%t.SUB.'0'))+ ANG
5634	datan2dt%t=datan2dt%t*RAD_TO_DEG_
5635	endif
5636	else !
5637	if(knob) then !knob 2
5638	call ass(datan2dt)
5639	call varfk2(S2)
5640	if(s1%kind==1) then !!
5641	ANG=ATAN2(varf2.SUB.'0',S1%R)
5642	if(pil>abs(ANG).or.abs(ANG)>pim) then !!!
5643	w%r=varf2/s1%r
5644	w=atan(w)
5645	datan2dt%t=w%r
5646	datan2dt%t=datan2dt%t-(datan2dt%t.SUB.'0')+ ANG
5647	datan2dt%t=datan2dt%t*RAD_TO_DEG_
5648	else !!!
5649	if((varf2.sub.'0')<0.0_dp) si=-1.0_dp
5650	w%r=s1%t/SQRT(varf22+s1%t2)
5651	w=acos(w)
5652	datan2dt%t=w%r
5653	datan2dt%t=si*(datan2dt%t-(datan2dt%t.SUB.'0'))+ ANG
5654	datan2dt%t=datan2dt%t*RAD_TO_DEG_
5655	endif !!!
5656	else !!
5657	ANG=ATAN2(varf2.SUB.'0',S1%t.sub.'0')
5658	if(pil>abs(ANG).or.abs(ANG)>pim) then !!!
5659	w%r=varf2/s1%t
5660	w=atan(w)
5661	datan2dt%t=w%r
5662	datan2dt%t=datan2dt%t-(datan2dt%t.SUB.'0')+ ANG
5663	datan2dt%t=datan2dt%t*RAD_TO_DEG_
5664	else !!!
5665	if((varf2.sub.'0')<0.0_dp) si=-1.0_dp
5666	w%r=s1%t/SQRT(varf22+s1%t2)
5667	w=acos(w)
5668	datan2dt%t=w%r
5669	datan2dt%t=si*(datan2dt%t-(datan2dt%t.SUB.'0'))+ ANG
5670	datan2dt%t=datan2dt%t*RAD_TO_DEG_
5671	endif !!!
5672	endif !!
5673
5674	else !knob 2
5675
5676	if(s1%kind==1) then !!
5677	datan2dt%R=ATAN2(S2%R,S1%R)*RAD_TO_DEG_
5678	datan2dt%kind=1
5679	else !!
5680	call ass(datan2dt)
5681	ANG=ATAN2(S2%R,S1%t.sub.'0')
5682	if(pil>abs(ANG).or.abs(ANG)>pim) then !!!
5683	w%r=S2%R/s1%t
5684	w=atan(w)
5685	datan2dt%t=w%r
5686	datan2dt%t=datan2dt%t-(datan2dt%t.SUB.'0')+ ANG
5687	datan2dt%t=datan2dt%t*RAD_TO_DEG_
5688	else !!!
5689	if((S2%R)<0.0_dp) si=-1.0_dp
5690	w%r=s1%t/SQRT(S2%R2+s1%t2)
5691	w=acos(w)
5692	datan2dt%t=w%r
5693	datan2dt%t=si*(datan2dt%t-(datan2dt%t.SUB.'0'))+ ANG
5694	datan2dt%t=datan2dt%t*RAD_TO_DEG_
5695	endif !!!
5696	endif !!
5697
5698	endif !knob 2
5699
5700
5701	endif !
5702	call kill(w)
5703	master=localmaster
5704	ELSE !1
5705	localmaster=master
5706	! call ass(datan2dt)
5707	call alloc(w)
5708	if(s1%kind==2) then !
5709	call ass(datan2dt)
5710	ANG=ATAN2(S2%R,S1%T.SUB.'0')
5711	if(pil>abs(ANG).or.abs(ANG)>pim) then
5712	w%r=s2%r/s1%t
5713	w=atan(w)
5714	datan2dt%t=w%r
5715	datan2dt%t=datan2dt%t-(datan2dt%t.SUB.'0')+ ANG
5716	datan2dt%t=datan2dt%t*RAD_TO_DEG_
5717	else
5718	if(s2%r<0.0_dp) si=-1.0_dp
5719	w%r=s1%t/SQRT(s2%r2+s1%t2)
5720	w=acos(w)
5721	datan2dt%t=w%r
5722	datan2dt%t=si*(datan2dt%t-(datan2dt%t.SUB.'0'))+ ANG
5723	datan2dt%t=datan2dt%t*RAD_TO_DEG_
5724	endif
5725	else !
5726	if(knob) then ! knob 3
5727	call ass(datan2dt)
5728	if(s2%kind==1) then !!
5729	call varfk1(S1)
5730	ANG=ATAN2(S2%R,varf1.SUB.'0')
5731	if(pil>abs(ANG).or.abs(ANG)>pim) then
5732	w%r=s2%r/varf1
5733	w=atan(w)
5734	datan2dt%t=w%r
5735	datan2dt%t=datan2dt%t-(datan2dt%t.SUB.'0')+ ANG
5736	datan2dt%t=datan2dt%t*RAD_TO_DEG_
5737	else
5738	if(s2%r<0.0_dp) si=-1.0_dp
5739	w%r=varf1/SQRT(s2%r2+varf12)
5740	w=acos(w)
5741	datan2dt%t=w%r
5742	datan2dt%t=si*(datan2dt%t-(datan2dt%t.SUB.'0'))+ ANG
5743	datan2dt%t=datan2dt%t*RAD_TO_DEG_
5744	endif
5745	else !! s2%kind ==2
5746	call varfk1(S1)
5747	ANG=ATAN2(S2%t.sub.'0',varf1.SUB.'0')
5748	if(pil>abs(ANG).or.abs(ANG)>pim) then
5749	w%r=s2%T/varf1
5750	w=atan(w)
5751	datan2dt%t=w%r
5752	datan2dt%t=datan2dt%t-(datan2dt%t.SUB.'0')+ ANG
5753	datan2dt%t=datan2dt%t*RAD_TO_DEG_
5754	else
5755	if((S2%t.sub.'0')<0.0_dp) si=-1.0_dp
5756	w%r=varf1/SQRT(s2%T2+varf12)
5757	w=acos(w)
5758	datan2dt%t=w%r
5759	datan2dt%t=si*(datan2dt%t-(datan2dt%t.SUB.'0'))+ ANG
5760	datan2dt%t=datan2dt%t*RAD_TO_DEG_
5761	endif
5762	endif !!
5763	else ! knob 3
5764	if(s2%kind==1) then !!
5765	datan2dt%R=ATAN2(S2%R,S1%R)*RAD_TO_DEG_
5766	datan2dt%kind=1
5767	else !! s2%kind ==2
5768	call ass(datan2dt)
5769	ANG=ATAN2(S2%t.sub.'0',S1%r)
5770	if(pil>abs(ANG).or.abs(ANG)>pim) then
5771	w%r=s2%T/S1%r
5772	w=atan(w)
5773	datan2dt%t=w%r
5774	datan2dt%t=datan2dt%t-(datan2dt%t.SUB.'0')+ ANG
5775	datan2dt%t=datan2dt%t*RAD_TO_DEG_
5776	else
5777	if((S2%t.sub.'0')<0.0_dp) si=-1.0_dp
5778	w%r=S1%r/SQRT(s2%T2+S1%r2)
5779	w=acos(w)
5780	datan2dt%t=w%r
5781	datan2dt%t=si*(datan2dt%t-(datan2dt%t.SUB.'0'))+ ANG
5782	datan2dt%t=datan2dt%t*RAD_TO_DEG_
5783	endif
5784	endif !!
5785	endif ! knob 3
5786	endif !
5787
5788	call kill(w)
5789	master=localmaster
5790	endif
5791
5792	END FUNCTION datan2dt
5793
5794	FUNCTION absoftdatan2dr( S2,S1 )
5795	implicit none
5796	real(dp) absoftdatan2dr
5797	real(dp), INTENT (IN) :: S2,S1
5798	absoftdatan2dr=ATAN2(S2,S1)*RAD_TO_DEG_
5799	END FUNCTION absoftdatan2dr
5800
5801	FUNCTION dasint( S1 )
5802	implicit none
5803	TYPE (real_8) dasint
5804	TYPE (real_8), INTENT (IN) :: S1
5805	TYPE (taylor) w
5806	integer localmaster
5807
5808	select case(s1%kind)
5809	case(m1)
5810	dasint%r=arcsin(s1%r)
5811	! dasint%r=asin(s1%r)
5812	dasint%kind=1
5813	case(m2)
5814	localmaster=master
5815	call ass(dasint)
5816	call alloc(w)
5817	w=s1%t
5818	w= asin(w)
5819	dasint%t=w
5820	call kill(w)
5821	master=localmaster
5822	case(m3)
5823	if(knob) then
5824	localmaster=master
5825	call ass(dasint)
5826	call alloc(w)
5827	call varfk1(S1)
5828	w=varf1
5829	w= asin(w)
5830	dasint%t=w
5831	call kill(w)
5832	master=localmaster
5833	else
5834	dasint%r=asin(s1%r)
5835	dasint%kind=1
5836	endif
5837	! call unvarkind3(S1)
5838	end select
5839	END FUNCTION dasint
5840
5841	FUNCTION dacost( S1 )
5842	implicit none
5843	TYPE (real_8) dacost
5844	TYPE (real_8), INTENT (IN) :: S1
5845	TYPE (taylor) w
5846	integer localmaster
5847
5848	select case(s1%kind)
5849	case(m1)
5850	dacost%r=arccos(s1%r)
5851	dacost%kind=1
5852	case(m2)
5853	localmaster=master
5854	call ass(dacost)
5855	call alloc(w)
5856	w=s1%t
5857	w= acos(w)
5858	dacost%t=w
5859	call kill(w)
5860	master=localmaster
5861	case(m3)
5862	if(knob) then
5863	localmaster=master
5864	call ass(dacost)
5865	call alloc(w)
5866	call varfk1(S1)
5867	w=varf1
5868	w= acos(w)
5869	dacost%t=w
5870	call kill(w)
5871	master=localmaster
5872	else
5873	dacost%r=acos(s1%r)
5874	dacost%kind=1
5875	endif
5876
5877	case default
5878	w_p=0
5879	w_p%nc=2
5880	w_p%fc='((1X,A72,/,1x,a72))'
5881	w_p%fi='(2((1X,i4)))'
5882	w_p%c(1)= " trouble in dacost "
5883	w_p%c(2)= "s1%kind "
5884	w_p=(/s1%kind /)
5885	! call !write_e(0)
5886	end select
5887	END FUNCTION dacost
5888
5889
5890	FUNCTION dcosht( S1 )
5891	implicit none
5892	TYPE (real_8) dcosht
5893	TYPE (real_8), INTENT (IN) :: S1
5894	TYPE (complextaylor) w
5895	integer localmaster
5896
5897	select case(s1%kind)
5898	case(m1)
5899	dcosht%r=cosh(s1%r)
5900	dcosht%kind=1
5901	case(m2)
5902	localmaster=master
5903	call ass(dcosht)
5904	call alloc(w)
5905	w%r=s1%t
5906	w= cosh(w)
5907	dcosht%t=w%r
5908	call kill(w)
5909	master=localmaster
5910	case(m3)
5911	if(knob) then
5912	localmaster=master
5913	call ass(dcosht)
5914	call alloc(w)
5915	call varfk1(S1)
5916	w%r=varf1
5917	w= cosh(w)
5918	dcosht%t=w%r
5919	call kill(w)
5920	master=localmaster
5921	else
5922	dcosht%r=cosh(s1%r)
5923	dcosht%kind=1
5924	endif
5925	case default
5926	w_p=0
5927	w_p%nc=2
5928	w_p%fc='((1X,A72,/,1x,a72))'
5929	w_p%fi='(2((1X,i4)))'
5930	w_p%c(1)= " trouble in dcosht "
5931	w_p%c(2)= "s1%kind "
5932	w_p=(/s1%kind /)
5933	! call !write_e(0)
5934	end select
5935	END FUNCTION dcosht
5936
5937	FUNCTION dsinht( S1 )
5938	implicit none
5939	TYPE (real_8) dsinht
5940	TYPE (real_8), INTENT (IN) :: S1
5941	TYPE (complextaylor) w
5942	integer localmaster
5943
5944	select case(s1%kind)
5945	case(m1)
5946	dsinht%r=sinh(s1%r)
5947	dsinht%kind=1
5948	case(m2)
5949	localmaster=master
5950	call ass(dsinht)
5951	call alloc(w)
5952	w%r=s1%t
5953	w= sinh(w)
5954	dsinht%t=w%r
5955	call kill(w)
5956	master=localmaster
5957	case(m3)
5958	if(knob) then
5959	localmaster=master
5960	call ass(dsinht)
5961	call alloc(w)
5962	call varfk1(S1)
5963	w%r=varf1
5964	w= sinh(w)
5965	dsinht%t=w%r
5966	call kill(w)
5967	master=localmaster
5968	else
5969	dsinht%r=sinh(s1%r)
5970	dsinht%kind=1
5971	endif
5972
5973	case default
5974	w_p=0
5975	w_p%nc=2
5976	w_p%fc='((1X,A72,/,1x,a72))'
5977	w_p%fi='(2((1X,i4)))'
5978	w_p%c(1)= " trouble in dsinht "
5979	w_p%c(2)= "s1%kind "
5980	w_p=(/s1%kind /)
5981	! call !write_e(0)
5982	end select
5983	END FUNCTION dsinht
5984
5985	FUNCTION dtanht( S1 )
5986	implicit none
5987	TYPE (real_8) dtanht
5988	TYPE (real_8), INTENT (IN) :: S1
5989	integer localmaster
5990
5991	select case(s1%kind)
5992	case(m1)
5993	dtanht%r=tanh(s1%r)
5994	dtanht%kind=1
5995	case(m2)
5996	localmaster=master
5997	call ass(dtanht)
5998	dtanht%t=tanh(s1%t)
5999	master=localmaster
6000	case(m3)
6001	if(knob) then
6002	localmaster=master
6003	call ass(dtanht)
6004	call varfk1(S1)
6005	dtanht%t=tanh(varf1)
6006	master=localmaster
6007	else
6008	dtanht%r=tanh(s1%r)
6009	dtanht%kind=1
6010	endif
6011
6012	case default
6013	w_p=0
6014	w_p%nc=2
6015	w_p%fc='((1X,A72,/,1x,a72))'
6016	w_p%fi='(2((1X,i4)))'
6017	w_p%c(1)= " trouble in dtanht "
6018	w_p%c(2)= "s1%kind "
6019	w_p=(/s1%kind /)
6020	! call !write_e(0)
6021	end select
6022	END FUNCTION dtanht
6023
6024	SUBROUTINE MAPreal_8(S2,S1)
6025	implicit none
6026	type (damap),INTENT(inout)::S2
6027	type (real_8),INTENT(IN),dimension(:)::S1
6028	integer i
6029
6030	call check_snake
6031
6032	do i=1,nd2
6033	s2%v(i)=s1(i) !%t
6034	enddo
6035	END SUBROUTINE MAPreal_8
6036
6037	SUBROUTINE normal_p(S2,S1)
6038	implicit none
6039	type (normalform),INTENT(inOUT)::S2
6040	type (real_8),INTENT(IN),dimension(:)::S1
6041	type (damap) t
6042	integer i
6043	call alloc(t)
6044
6045	call check_snake
6046
6047	do i=1,nd2
6048	t%v(i)=s1(i) !%t
6049	enddo
6050	s2=t
6051
6052	call kill(t)
6053	END SUBROUTINE normal_p
6054
6055
6056	SUBROUTINE real_8MAP(S1,S2)
6057	implicit none
6058	type (damap),INTENT(in)::S2
6059	type (real_8),INTENT(inout),dimension(:)::S1
6060	integer i
6061
6062	call check_snake
6063
6064	do i=1,nd2
6065	s1(i)=s2%v(i)
6066	enddo
6067
6068	END SUBROUTINE real_8MAP
6069
6070
6071	! radiation
6072	SUBROUTINE ENVreal_8(S2,S1)
6073	implicit none
6074	type (ENV_8),INTENT(inout),dimension(:)::S2
6075	type (real_8),INTENT(IN),dimension(:)::S1
6076	integer i
6077
6078	call check_snake
6079
6080	do i=1,nd2
6081	s2(i)%v=s1(i) !%t
6082	enddo
6083	END SUBROUTINE ENVreal_8
6084
6085	SUBROUTINE real_8ENV(S1,S2)
6086	implicit none
6087	type (ENV_8),INTENT(in),dimension(:)::S2
6088	type (real_8),INTENT(inout),dimension(:)::S1
6089	integer i
6090
6091	call check_snake
6092
6093	do i=1,nd2
6094	s1(i)=s2(i)%v
6095	enddo
6096	END SUBROUTINE real_8ENV
6097
6098	SUBROUTINE ENV_8RADTAYLOR(S2,S1)
6099	implicit none
6100	type (ENV_8),INTENT(inout),dimension(:)::S2
6101	type (RADTAYLOR),INTENT(IN),dimension(:)::S1
6102	integer i,J
6103
6104	call check_snake
6105
6106	do i=1,nd2
6107	s2(i)%v=s1(i)%V !%t
6108	enddo
6109	do i=1,nd2
6110	do J=1,nd2
6111	s2(i)%E(J)=s1(i)%E(J) !%t
6112	enddo
6113	enddo
6114
6115	END SUBROUTINE ENV_8RADTAYLOR
6116
6117	SUBROUTINE RADTAYLORENV_8(S1,S2)
6118	implicit none
6119	type (ENV_8),INTENT(in),dimension(:)::S2
6120	type (RADTAYLOR),INTENT(inout),dimension(:)::S1
6121	integer i,J
6122
6123	call check_snake
6124
6125	do i=1,nd2
6126	s1(i)%V= s2(i)%v !%t
6127	enddo
6128	do i=1,nd2
6129	do J=1,nd2
6130	s1(i)%E(J)=s2(i)%E(J) !%t
6131	enddo
6132	enddo
6133
6134	END SUBROUTINE RADTAYLORENV_8
6135
6136	SUBROUTINE RADTAYLORprobe_8(S1,S2)
6137	implicit none
6138	type (probe_8),INTENT(in) ::S2
6139	type (RADTAYLOR),INTENT(inout),dimension(:)::S1
6140	integer i,J
6141
6142	call check_snake
6143
6144	do i=1,nd2
6145	s1(i)%V= s2%X(I) !%t
6146	enddo
6147	do i=1,nd2
6148	do J=1,nd2
6149	s1(i)%E(J)=s2%E_IJ(i,J) !%t
6150	enddo
6151	enddo
6152
6153	END SUBROUTINE RADTAYLORprobe_8
6154
6155	SUBROUTINE beamprobe_8(S2,Sprobe_8)
6156	implicit none
6157	type (beamenvelope),INTENT(inout)::S2
6158	type (radtaylor) S1(ndim2)
6159	type (probe_8),INTENT(IN)::Sprobe_8
6160
6161	call alloc(s1,nd2)
6162
6163	call check_snake
6164
6165	S1=Sprobe_8
6166
6167
6168	s2=s1
6169
6170	call kill(s1,nd2)
6171
6172
6173	END SUBROUTINE beamprobe_8
6174
6175	SUBROUTINE resetenv(S2)
6176	implicit none
6177	type (env_8),INTENT(INOUT)::S2
6178	integer i
6179
6180	call resetpoly(s2%v)
6181	do i=1,nd2
6182	call resetpoly(s2%e(i))
6183	call resetpoly(s2%sigma0(i))
6184	call resetpoly(s2%sigmaf(i))
6185	enddo
6186
6187	END SUBROUTINE resetenv
6188
6189
6190	SUBROUTINE allocenv(S2)
6191	implicit none
6192	type (env_8),INTENT(INOUT)::S2
6193	integer i
6194
6195	call alloc(s2%v)
6196	do i=1,nd2
6197	call alloc(s2%e(i))
6198	call alloc(s2%sigma0(i))
6199	call alloc(s2%sigmaf(i))
6200	enddo
6201
6202	END SUBROUTINE allocenv
6203
6204	SUBROUTINE allocenvn(S2,K)
6205	implicit none
6206	type (env_8),INTENT(INOUT),dimension(:)::S2
6207	INTEGER,optional,INTENT(IN)::k
6208	INTEGER J,i,N
6209
6210	if(present(k)) then
6211	I=LBOUND(S2,DIM=1)
6212	N=LBOUND(S2,DIM=1)+K-1
6213	else
6214	I=LBOUND(S2,DIM=1)
6215	N=UBOUND(S2,DIM=1)
6216	endif
6217
6218	DO J=I,N
6219	call allocenv(s2(j))
6220	enddo
6221
6222
6223	END SUBROUTINE allocenvn
6224
6225	SUBROUTINE e_OPT(S1,S2,s3,s4,s5,s6,s7,s8,s9,s10)
6226	implicit none
6227	type (env_8),INTENT(INout)::S1
6228	type (env_8),optional, INTENT(INout):: S2,s3,s4,s5,s6,s7,s8,s9,s10
6229	call allocenv(s1)
6230	if(present(s2)) call allocenv(s2)
6231	if(present(s3)) call allocenv(s3)
6232	if(present(s4)) call allocenv(s4)
6233	if(present(s5)) call allocenv(s5)
6234	if(present(s6)) call allocenv(s6)
6235	if(present(s7)) call allocenv(s7)
6236	if(present(s8)) call allocenv(s8)
6237	if(present(s9)) call allocenv(s9)
6238	if(present(s10))call allocenv(s10)
6239	END SUBROUTINE e_opt
6240
6241
6242	! SUBROUTINE allocenvn(S2,n)
6243	! implicit none
6244	! type (env_8),INTENT(INOUT),dimension(:)::S2
6245	! integer, INTENT(IN)::n
6246	! integer i
6247	!
6248	! do i=1,n
6249	! call allocenv(s2(i))
6250	! enddo
6251
6252	! END SUBROUTINE allocenvn
6253
6254	SUBROUTINE resetenvn(S2,K)
6255	implicit none
6256	type (env_8),INTENT(INOUT),dimension(:)::S2
6257	INTEGER,optional,INTENT(IN)::k
6258	INTEGER J,i,N
6259
6260	if(present(k)) then
6261	I=LBOUND(S2,DIM=1)
6262	N=LBOUND(S2,DIM=1)+K-1
6263	else
6264	I=LBOUND(S2,DIM=1)
6265	N=UBOUND(S2,DIM=1)
6266	endif
6267
6268	DO J=I,N
6269	call resetenv(s2(j))
6270	enddo
6271
6272
6273	END SUBROUTINE resetenvn
6274
6275	! SUBROUTINE resetenvn(S2,n)
6276	! implicit none
6277	! type (env_8),INTENT(INOUT),dimension(:)::S2
6278	! integer, INTENT(IN)::n
6279	! integer i
6280	!
6281	! do i=1,n
6282	! call resetenv(s2(i))
6283	! enddo
6284	!
6285	!
6286	! END SUBROUTINE resetenvn
6287
6288
6289
6290	SUBROUTINE beamENV_8(S2,SENV_8)
6291	implicit none
6292	type (beamenvelope),INTENT(inout)::S2
6293	type (radtaylor) S1(ndim2)
6294	! type (ENV_8),INTENT(IN)::SENV_8(ndim2)
6295	type (ENV_8),INTENT(IN)::SENV_8(6)
6296
6297	call alloc(s1,nd2)
6298
6299	call check_snake
6300
6301	S1=SENV_8
6302
6303
6304	s2=s1
6305
6306	call kill(s1,nd2)
6307
6308
6309	END SUBROUTINE beamENV_8
6310
6311
6312	SUBROUTINE varfk1(S2)
6313	implicit none
6314	! type (real_8),INTENT(INOUT)::S2
6315	type (real_8),INTENT(IN):: S2
6316
6317	if(knob) then
6318	if(nb_==0) then
6319	varf1=(/S2%R,S2%S/).var.(s2%i+npara_fpp)
6320	elseif(s2%nb==nb_) then
6321	varf1=(/S2%R,S2%S/).var.(s2%i+npara_fpp-s2%g+1)
6322	else
6323	varf1=S2%R
6324	endif
6325	else ! Not a knob
6326	stop 333
6327	varf1=(/S2%R,0.0_dp/).var.0 ! this is a buggy line never used
6328	endif
6329
6330
6331	end SUBROUTINE varfk1
6332
6333	SUBROUTINE varfk2(S2)
6334	implicit none
6335	! type (real_8),INTENT(INOUT)::S2
6336	type (real_8),INTENT(IN):: S2
6337
6338	if(knob) then
6339	if(nb_==0) then
6340	varf2=(/S2%R,S2%S/).var.(s2%i+npara_fpp)
6341	elseif(s2%nb==nb_) then
6342	varf2=(/S2%R,S2%S/).var.(s2%i+npara_fpp-s2%g+1)
6343	else
6344	varf2=S2%R
6345	endif
6346	else ! Not a knob
6347	stop 334
6348	varf2=(/S2%R,0.0_dp/).var.0 ! this is a buggy line never used
6349	endif
6350
6351
6352	end SUBROUTINE varfk2
6353
6354
6355	! SINH(X)/X DONE partly COSY-INFINITY WISE FOR TPSA
6356	function SINH_HR(X)
6357	implicit none
6358	real(dp), INTENT (IN) :: X
6359	real(dp) SINH_HR,Y,NORM0,NORM,SINH_HR0
6360	logical(lp) NOTDONE,CHECK
6361	INTEGER I,ipause,mypauses
6362
6363	if(abs(x)<sinhx_x_min) then
6364	NOTDONE=.TRUE.
6365	CHECK=.TRUE.
6366
6367	SINH_HR=1.0_dp
6368	Y=1.0_dp
6369	I=1
6370	NORM0=1e5_dp
6371	DO WHILE(I<NMAX_pol.AND.NOTDONE)
6372	Y=YX*2/REAL(I+1,kind=DP)/REAL(I+2,kind=DP)
6373	SINH_HR0=SINH_HR+Y
6374	NORM=ABS(SINH_HR-SINH_HR0)
6375	IF(NORM<=EPS_real_poly.AND.CHECK) THEN
6376	NORM0=NORM
6377	CHECK=.FALSE.
6378	ELSE
6379	IF(NORM>=NORM0) THEN
6380	NOTDONE=.FALSE.
6381	ELSE
6382	NORM0=NORM
6383	ENDIF
6384	ENDIF
6385	SINH_HR=SINH_HR0
6386	I=I+2
6387	ENDDO
6388	IF(I==NMAX_pol) THEN
6389	line="NO CONVERGENCE IN SINH_HR"
6390	ipause=mypauses(NMAX_pol,line)
6391	ENDIF
6392	else
6393	SINH_HR=SINH(X)/X
6394	endif
6395	return
6396	end function SINH_HR
6397
6398
6399	! SIN(X)/X DONE partly COSY-INFINITY WISE FOR TPSA
6400	function SIN_HR(X)
6401	implicit none
6402	real(dp), INTENT (IN) :: X
6403	real(dp) SIN_HR,Y,NORM0,NORM,SINH_HR0
6404	logical(lp) NOTDONE,CHECK
6405	INTEGER I,ipause,mypauses
6406
6407	if(abs(x)<sinhx_x_min) then
6408	NOTDONE=.TRUE.
6409	CHECK=.TRUE.
6410
6411	SIN_HR=1.0_dp
6412	Y=1.0_dp
6413	I=1
6414	NORM0=1e5_dp
6415	DO WHILE(I<NMAX_pol.AND.NOTDONE)
6416	Y=-YX*2/REAL(I+1,kind=DP)/REAL(I+2,kind=DP)
6417	SINH_HR0=SIN_HR+Y
6418	NORM=ABS(SIN_HR-SINH_HR0)
6419	IF(NORM<=EPS_real_poly.AND.CHECK) THEN
6420	NORM0=NORM
6421	CHECK=.FALSE.
6422	ELSE
6423	IF(NORM>=NORM0) THEN
6424	NOTDONE=.FALSE.
6425	ELSE
6426	NORM0=NORM
6427	ENDIF
6428	ENDIF
6429	SIN_HR=SINH_HR0
6430	I=I+2
6431	ENDDO
6432	IF(I==NMAX_pol) THEN
6433	line="NO CONVERGENCE IN SINH_HR"
6434	ipause=mypauses(NMAX_pol,line)
6435	ENDIF
6436	else
6437	SIN_HR=SIN(X)/X
6438	endif
6439	return
6440	end function SIN_HR
6441
6442	FUNCTION sinX_Xt( S1 )
6443	implicit none
6444	integer ipause, mypauses
6445	TYPE (real_8) sinX_Xt
6446	TYPE (real_8), INTENT (IN) :: S1
6447	integer localmaster
6448	TYPE(TAYLOR)SIN_HP
6449	TYPE(TAYLOR)Y,SINH_HP0
6450	real(dp) NORM0,NORM
6451	logical(lp) NOTDONE,CHECK
6452	INTEGER I
6453
6454	select case(s1%kind)
6455	case(m1)
6456	sinX_Xt%r=SIN_HR(s1%r)
6457	sinX_Xt%kind=1
6458	case(m2)
6459	CALL ALLOC(SIN_HP); CALL ALLOC(Y); CALL ALLOC(SINH_HP0);
6460	localmaster=master
6461	call ass(sinX_Xt)
6462
6463	if(ABS(S1%T.SUB.'0')<sinhx_x_minp) then
6464	NOTDONE=.TRUE.
6465	CHECK=.TRUE.
6466
6467	SIN_HP=1.0_dp
6468	Y=1.0_dp
6469	I=1
6470	NORM0=1e5_dp
6471	DO WHILE(I<NMAX_pol.AND.NOTDONE)
6472	Y=-YS1%T*2/REAL(I+1,kind=DP)/REAL(I+2,kind=DP)
6473	SINH_HP0=SIN_HP+Y
6474	NORM=full_abs(SIN_HP-SINH_HP0)
6475	IF(NORM<=EPS_real_poly.AND.CHECK) THEN
6476	NORM0=NORM
6477	CHECK=.FALSE.
6478	ELSE
6479	IF(NORM>=NORM0.and.(.not.check)) THEN ! sagan
6480	NOTDONE=.FALSE.
6481	ELSE
6482	NORM0=NORM
6483	ENDIF
6484	ENDIF
6485	SIN_HP=SINH_HP0
6486	I=I+2
6487	ENDDO
6488	IF(I==NMAX_pol) THEN
6489	line="NO CONVERGENCE IN SIN_HP"
6490	ipause=mypauses(NMAX_pol,line)
6491	ENDIF
6492	else
6493	SIN_HP=SIN(S1%T)/S1%T
6494	endif
6495	sinX_Xt%t= SIN_HP
6496
6497	master=localmaster
6498	CALL KILL(SIN_HP); CALL KILL(Y); CALL KILL(SINH_HP0);
6499	case(m3)
6500	if(knob) then
6501	localmaster=master
6502	call ass(sinX_Xt)
6503	call varfk1(S1)
6504	if(ABS(varf1.SUB.'0')<sinhx_x_minp) then
6505	NOTDONE=.TRUE.
6506	CHECK=.TRUE.
6507
6508	SIN_HP=1.0_dp
6509	Y=1.0_dp
6510	I=1
6511	NORM0=1e5_dp
6512	DO WHILE(I<NMAX_pol.AND.NOTDONE)
6513	Y=-Yvarf1*2/REAL(I+1,kind=DP)/REAL(I+2,kind=DP)
6514	SINH_HP0=SIN_HP+Y
6515	NORM=full_abs(SIN_HP-SINH_HP0)
6516	IF(NORM<=EPS_real_poly.AND.CHECK) THEN
6517	NORM0=NORM
6518	CHECK=.FALSE.
6519	ELSE
6520	IF(NORM>=NORM0.and.(.not.check)) THEN ! sagan
6521	NOTDONE=.FALSE.
6522	ELSE
6523	NORM0=NORM
6524	ENDIF
6525	ENDIF
6526	SIN_HP=SINH_HP0
6527	I=I+2
6528	ENDDO
6529	IF(I==NMAX_pol) THEN
6530	line="NO CONVERGENCE IN SINH_HP"
6531	ipause=mypauses(NMAX_pol,line)
6532	ENDIF
6533	else
6534	SIN_HP=SIN(varf1)/varf1
6535	endif
6536	sinX_Xt%t= SIN_HP
6537	master=localmaster
6538	else
6539	sinX_Xt%r= SIN_HR(S1%r)
6540	sinX_Xt%kind=1
6541	endif
6542	case default
6543	w_p=0
6544	w_p%nc=2
6545	w_p%fc='((1X,A72,/,1x,a72))'
6546	w_p%fi='(2((1X,i4)))'
6547	w_p%c(1)= " trouble in sinX_XT "
6548	w_p%c(2)= "s1%kind "
6549	w_p=(/s1%kind /)
6550	! call !write_e(0)
6551	end select
6552	END FUNCTION sinX_XT
6553
6554	FUNCTION sinHX_Xt( S1 )
6555	implicit none
6556	integer ipause, mypauses
6557	TYPE (real_8) sinHX_Xt
6558	TYPE (real_8), INTENT (IN) :: S1
6559	integer localmaster
6560	TYPE(TAYLOR)SIN_HP
6561	TYPE(TAYLOR)Y,SINH_HP0
6562	real(dp) NORM0,NORM
6563	logical(lp) NOTDONE,CHECK
6564	INTEGER I
6565
6566	select case(s1%kind)
6567	case(m1)
6568	sinHX_Xt%r=SINH_HR(s1%r)
6569	sinHX_Xt%kind=1
6570	case(m2)
6571	CALL ALLOC(SIN_HP); CALL ALLOC(Y); CALL ALLOC(SINH_HP0);
6572	localmaster=master
6573	call ass(sinHX_Xt)
6574
6575	if(ABS(S1%T.SUB.'0')<sinhx_x_minp) then
6576	NOTDONE=.TRUE.
6577	CHECK=.TRUE.
6578
6579	SIN_HP=1.0_dp
6580	Y=1.0_dp
6581	I=1
6582	NORM0=1e5_dp
6583	DO WHILE(I<NMAX_pol.AND.NOTDONE)
6584	Y=YS1%T*2/REAL(I+1,kind=DP)/REAL(I+2,kind=DP)
6585	SINH_HP0=SIN_HP+Y
6586	NORM=full_abs(SIN_HP-SINH_HP0)
6587	IF(NORM<=EPS_real_poly.AND.CHECK) THEN
6588	NORM0=NORM
6589	CHECK=.FALSE.
6590	ELSE
6591	IF(NORM>=NORM0.and.(.not.check)) THEN ! sagan
6592	NOTDONE=.FALSE.
6593	ELSE
6594	NORM0=NORM
6595	ENDIF
6596	ENDIF
6597	SIN_HP=SINH_HP0
6598	I=I+2
6599	ENDDO
6600	IF(I==NMAX_pol) THEN
6601	line="NO CONVERGENCE IN SINH_HP"
6602	ipause=mypauses(NMAX_pol,line)
6603	ENDIF
6604	else
6605
6606	SIN_HP=SINH(S1%T)/S1%T
6607	endif
6608	sinHX_Xt%t= SIN_HP
6609
6610	master=localmaster
6611	CALL KILL(SIN_HP); CALL KILL(Y); CALL KILL(SINH_HP0);
6612	case(m3)
6613	if(knob) then
6614	localmaster=master
6615	call ass(sinHX_Xt)
6616	call varfk1(S1)
6617	if(ABS(varf1.SUB.'0')<sinhx_x_minp) then
6618	NOTDONE=.TRUE.
6619	CHECK=.TRUE.
6620
6621	SIN_HP=1.0_dp
6622	Y=1.0_dp
6623	I=1
6624	NORM0=1e5_dp
6625	DO WHILE(I<NMAX_pol.AND.NOTDONE)
6626	Y=Yvarf1*2/REAL(I+1,kind=DP)/REAL(I+2,kind=DP)
6627	SINH_HP0=SIN_HP+Y
6628	NORM=full_abs(SIN_HP-SINH_HP0)
6629	IF(NORM<=EPS_real_poly.AND.CHECK) THEN
6630	NORM0=NORM
6631	CHECK=.FALSE.
6632	ELSE
6633	IF(NORM>=NORM0.and.(.not.check)) THEN ! sagan
6634	NOTDONE=.FALSE.
6635	ELSE
6636	NORM0=NORM
6637	ENDIF
6638	ENDIF
6639	SIN_HP=SINH_HP0
6640	I=I+2
6641	ENDDO
6642	IF(I==NMAX_pol) THEN
6643	line="NO CONVERGENCE IN SINH_HP"
6644	ipause=mypauses(NMAX_pol,line)
6645	ENDIF
6646	else
6647	SIN_HP=SINH(varf1)/varf1
6648	endif
6649	sinHX_Xt%t= SIN_HP
6650	master=localmaster
6651	else
6652	sinHX_Xt%r= SINH_HR(S1%r)
6653	sinHX_Xt%kind=1
6654	endif
6655	case default
6656	w_p=0
6657	w_p%nc=2
6658	w_p%fc='((1X,A72,/,1x,a72))'
6659	w_p%fi='(2((1X,i4)))'
6660	w_p%c(1)= " trouble in sinHX_Xt "
6661	w_p%c(2)= "s1%kind "
6662	w_p=(/s1%kind /)
6663	! call !write_e(0)
6664	end select
6665	END FUNCTION sinHX_Xt
6666
6667	! remove small numbers
6668
6669	SUBROUTINE clean_real_8(S1,S2,prec)
6670	implicit none
6671	type (real_8),INTENT(INOUT)::S2
6672	type (real_8), intent(INOUT):: s1
6673	real(dp) prec
6674	integer i
6675	type(real_8) t
6676
6677	call alloc(t)
6678	t=s1
6679
6680	select case(s1%kind)
6681	case(m1)
6682	if(abs(t%r)<prec) t%r=0.0_dp
6683	case(m2)
6684	call clean_taylor(t%t,t%t,prec)
6685	case(m3)
6686	Write(6,*) " cannot clean a knob "
6687	stop 601
6688	case default
6689	w_p=0
6690	w_p%nc=2
6691	w_p%fc='((1X,A72,/,1x,a72))'
6692	w_p%fi='(2((1X,i4)))'
6693	w_p%c(1)= " trouble in clean_real_8 "
6694	w_p%c(2)= "s1%kind "
6695	w_p=(/s1%kind /)
6696	! call !write_e(0)
6697	end select
6698	s2=t
6699	call kill(t)
6700
6701
6702	END SUBROUTINE clean_real_8
6703
6704	SUBROUTINE flip_real_8(S1,S2,i)
6705	implicit none
6706	type (real_8),INTENT(INOUT)::S2
6707	type (real_8), intent(INOUT):: s1
6708	integer i
6709	if(s1%kind==2) then
6710	s2=s1 ! to make s2 taylor!
6711	call flip_taylor(S1%t,S2%t,i)
6712	endif
6713	end SUBROUTINE flip_real_8
6714
6715
6716	end module polymorphic_taylor

Note: See TracBrowser for help on using the repository browser.

Context Navigation

source: PSPA/madxPSPA/libs/ptc/src/m_real_polymorph.f90 @ 430

Download in other formats: