-- development test used to develope GENERIC_COMPLEX_ELEMENTARY_FUNCTIONS -- this test uses identities as a "sanity check" on algorithms. -- this is NOT an accuracy check -- this test excludes specific values that fail when there is not a signed zero -- generally all conjugate tests fail on the branch cut without a signed zero -- written 4/20/90 by Jon Squire -- revised 5/31/91 by Jon Squire - add more tests -- revised 7/17/91 by Jon Squire - fine tune some tests -- revised 9/1/91 by Jon Squire - add constraint tests -- revised 5/15/92 by Jon Squire - add more tests -- revised 1/17/93 by Jon Squire - fix instantiations to new std -- revised 9/20/96 by Jon Squire - fix IO, be ready for Ada '95 -- The actual tests for each base type are created by the -- instantiations at the end. Ignore compilation errors for types -- that do not exists. -- The instantiated procedures are: -- TEST_COMPLEX_IDENTITIES -- TEST_LONG_COMPLEX_IDENTITIES with text_io; use text_io; with generic_real_arrays; with generic_complex_types; with generic_complex_arrays; with complex_io; with generic_complex_elementary_functions; with generic_elementary_functions; with generic_complex_constraints; generic type real is digits <>; procedure test_generic_complex_identities; procedure test_generic_complex_identities is package real_arrays is new generic_real_arrays(real); use real_arrays; -- make operators visible for user -- make subprograms visible for instantiations package complex_types is new generic_complex_types(real); use complex_types; -- make operators visible for user -- make subprograms visible for instantiations package cx_io is new complex_io(real, complex); use cx_io; package complex_elementary_functions is new generic_complex_elementary_functions(real, complex, imaginary); use complex_elementary_functions; package complex_arrays is new generic_complex_arrays (real, real_vector, real_matrix, complex); use complex_arrays; -- make operators visible for user -- make subprograms visible for instantiations package real_io is new float_io(real); use real_io; package complex_constraints is new generic_complex_constraints(real, complex); use complex_constraints; package real_functions is new generic_elementary_functions(real); package rf renames real_functions; argument_error : exception renames complex_types.argument_error; PI : constant := 3.14159_26535_89793_23846_26433_83279_50288_41971 ; Pi_2 : constant := Pi / 2.0; two_Pi : constant := 2.0 * Pi; K : array(integer(-44)..44) of real; P : COMPLEX_VECTOR ( 1..k'length*k'length); ii : integer; complex_zero : complex := (0.0,0.0); complex_one : complex := (1.0,0.0); complex_i : complex := (0.0,1.0); maxerr : real := 0.0; perr : complex := (0.0,0.0); myval : complex := (0.0,0.0); ident : complex := (0.0,0.0); base : complex := (0.0,0.0); iden : complex := (0.0,0.0); excp : integer := 0; excl : integer := 0; errc : integer := 0; z : complex; procedure track(point:complex; gcef_value: complex; ident_value: complex; print_all: boolean := false) is error : real; begin if modulus(gcef_value) /= 0.0 then error := modulus(ident_value - gcef_value)/modulus(gcef_value); else error := modulus(ident_value - gcef_value); end if; if error > maxerr then if print_all then put("Track at z= "); put(point); put(", err= "); put(error); new_line; end if; maxerr := error; perr := point; myval := gcef_value; ident := ident_value; end if; if error > 0.01 then errc := errc + 1; end if; end track; procedure print_track is begin put(" exception count "); put(integer'image(excp)); put(" excluded count "); put(integer'image(excl)); put(" error count "); put(integer'image(errc)); new_line; put(" maxerr "); put(maxerr,aft=>3); new_line; put(" z at maxerr"); put(perr); new_line; put(" base "); put(myval); new_line; put(" identity "); put(ident); new_line; new_line; maxerr := 0.0; excp := 0; excl := 0; errc := 0; perr := (0.0,0.0); myval := (0.0,0.0); ident := (0.0,0.0); end print_track; function square(z:complex) return complex is x : real; y : real; begin x := (RE(z)-IM(z)) * (RE(z)+IM(z)); y := RE(z)*IM(z); return compose_from_cartesian(x,y+y); end square; function csgn(x,y:real) return real is -- for z=x+iy 1.0 in right plane begin -- -1.0 in left palne, 0.0 at origin if x > 0.0 then return 1.0; elsif x < 0.0 then return -1.0; else -- x = 0.0 if y > 0.0 then return 1.0; elsif y < 0.0 then return -1.0; else -- y = 0.0 return 0.0; end if; end if; end csgn; function sss(x,y:real) return real is -- sqrt sum of squares begin -- sqrt(x^2 + y^2) if x = 0.0 and then y = 0.0 then return 0.0 ; end if; if abs x > abs y then return abs x * rf.sqrt( 1.0 + (y/x)*(y/x) ); else return abs y * rf.sqrt( 1.0 + (x/y)*(x/y) ); end if; end sss; function sds(x,y:real) return real is -- sqrt difference of squares begin -- sqrt(x^2 - y^2) if x = 0.0 and then y = 0.0 then return 0.0 ; elsif x = y then return 0.0; end if; if abs x > abs y then return abs x * rf.sqrt( 1.0 - (y/x)*(y/x) ); else return abs y * rf.sqrt( (x/y)*(x/y) - 1.0 ); end if; end sds; function local_arcsin(z:complex) return complex is x : real := re(z); y : real := im(z); q1 : real := sss(x+1.0,y)/2.0; q2 : real := sss(x-1.0,y)/2.0; begin if abs x < real'epsilon and then abs y < real'epsilon then return z; elsif abs x > 1.0/real'epsilon or else abs y > 1.0/real'epsilon then return compose_from_cartesian(0.0, csgn(y,-x)*rf.log(2.0*sss(x,y))); end if; return compose_from_cartesian(rf.arcsin(q1-q2), csgn(y,-x)*rf.log(q1+q2+sds(q1+q2,1.0))); end local_arcsin; function local_arccos(z:complex) return complex is x : real := re(z); y : real := im(z); q1 : real := sss(x+1.0,y)/2.0; q2 : real := sss(x-1.0,y)/2.0; begin if abs x < real'epsilon and then abs y < real'epsilon then return compose_from_cartesian(pi_2,0.0); elsif abs x > 1.0/real'epsilon or else abs y > 1.0/real'epsilon then return compose_from_cartesian(pi_2, csgn(-y,x)*rf.log(2.0*sss(x,y))); end if; return compose_from_cartesian(rf.arccos(q1-q2), csgn(-y,x)*rf.log(q1+q2+sds(q1+q2,1.0))); end local_arccos; begin -- The values below are used positive and negative in all combinations -- real and imaginary parts k(0) := 0.0; k(1) := 0.0000001; k(2) := 0.000001; k(3) := 0.00001; k(4) := 0.001; k(5) := 0.01; k(6) := 0.02; k(7) := 0.06; k(8) := 0.1; k(9) := 0.2; k(10) := 0.25; k(11) := 0.35; k(12) := 0.5; k(13) := 0.7071; k(14) := 0.785; k(15) := 0.786; k(16) := 0.97; k(17) := 0.999; k(18) := 1.0; k(19) := 1.001; k(20) := 1.1; k(21) := 1.01; k(22) := 1.25; k(23) := 1.35; k(24) := 1.41421; k(25) := 1.5; k(26) := 1.54; k(27) := 1.57079; k(28) := 1.63; k(29) := 2.0; k(30) := 2.2; k(31) := 2.4; k(32) := 2.71828; k(33) := 3.0; k(34) := 3.141; k(35) := 3.142; k(36) := 8.0; k(37) := 12.0; k(38) := 16.0; k(39) := 20.0; k(40) := 30.0; k(41) := 64.0; k(42) := 500.0; k(43) := 1000.0; k(44) := 2500.0; for i in 1..k'last loop k(-i) := -k(i); end loop; ii := p'first; for i in k'range loop for j in k'range loop P(ii) := compose_from_cartesian(k(i),k(j)); ii := ii+1; end loop; end loop; put_line("test complex identities"); put_line("test sqrt"); put_line("sqrt(z)=polar(sqrt(|z|),arg(z)/2)"); for i in p'range loop begin z := p(i); base:=sqrt(z); iden:=compose_from_polar(real_functions.sqrt(abs z), argument(z)/2.0); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("sqrt(z)**2=z"); for i in p'range loop begin z := p(i); base:=square(sqrt(z)); iden:=z; track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("sqrt(z**2)=z z in right half plane, not neg imag axis"); for i in p'range loop begin z := p(i); if RE(z) > 0.0 or else ( RE(z) = 0.0 and IM(z) >= 0.0 )then base := z; iden := sqrt(square(z)); track(z,base,iden); else excl := excl + 1; end if; exception when others => excp := excp + 1; end; end loop; print_track; put_line("sqrt(z)=sqrt(abs z)*exp(complex_i*argument(z)/2.0)"); for i in p'range loop begin z := p(i); base := sqrt(z); iden := real_functions.sqrt(abs z)*exp(complex_i*argument(z)/2.0); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("sqrt(1/z)=1/sqrt(z) not on negative real axis if no -0 "); for i in p'range loop begin z := p(i); if not ( IM(z)=0.0 and RE(z) < 0.0 ) then base := sqrt(1.0/z); iden := 1.0/sqrt(z); track(z,base,iden); else excl := excl + 1; end if; exception when others => excp := excp + 1; end; end loop; print_track; put_line("conjugate(sqrt(z)) = sqrt(conjugate(z)) not on branch cut"); for i in p'range loop begin z := p(i); if IM(z) /= 0.0 or RE(z) > 0.0 then base := conjugate(sqrt(z)); iden := sqrt(conjugate(z)); track(z,base,iden); else excl := excl + 1; end if; exception when others => excp := excp + 1; end; end loop; print_track; put_line("test exp"); put_line("exp(log(z))=z modulus z non zero"); for i in p'range loop begin z := p(i); if modulus(z) /= 0.0 then base := z; iden := exp(log(z)); track(z,base,iden); else excl := excl + 1; end if; exception when others => excp := excp + 1; end; end loop; print_track; put_line("exp(-z)=1/exp(z) not on negative real axis"); for i in p'range loop begin z := p(i); base := exp(-z); iden := 1.0/exp(z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("exp(z)=exp(z + i 2Pi) periodic"); for i in p'range loop begin z := p(i); base := exp(z); iden := exp(z+complex_i*two_pi); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("conjugate(exp(z)) = exp(conjugate(z))"); for i in p'range loop begin z := p(i); base := conjugate(exp(z)); iden := exp(conjugate(z)); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("test log"); put_line("log(exp(z))=z <- limit check z to im z in [-Pi,Pi]"); for i in p'range loop begin z := p(i); if IM(z)>-pi and IM(z) excp := excp + 1; end; end loop; print_track; put_line("log(1/z)=-log(z) not on negative real axis"); for i in p'range loop begin z := p(i); if not on_branch_cut_log(z) then base := log(1.0/z); iden := -log(z); track(z,base,iden); else excl := excl + 1; end if; exception when others => excp := excp + 1; end; end loop; print_track; put_line("conjugate(log(z)) = log(conjugate(z)) not on branch cut"); for i in p'range loop begin z := p(i); if not on_branch_cut_log(z) then base := conjugate(log(z)); iden := log(conjugate(z)); track(z,base,iden); else excl := excl + 1; end if; exception when others => excp := excp + 1; end; end loop; print_track; put_line("test sin"); put_line("sin(arcsin(z))=z"); for i in p'range loop begin z := p(i); base := z; iden := sin(arcsin(z)); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("sin(z)=(exp(iz)-exp(-iz))/2i"); for i in p'range loop begin z := p(i); base := sin(z); iden := (exp(complex_i*z)-exp(-complex_i*z))/(0.0,2.0); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("sin(z) = -i sinh iz"); for i in p'range loop begin z := p(i); base := sin(z); iden := -complex_i*sinh(complex_i*z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("sin(z) = cos(z - Pi/2)"); for i in p'range loop begin z := p(i); base := sin(z); iden := cos(z-Pi_2); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("sin(z) = -sin(-z)"); for i in p'range loop begin z := p(i); base := sin(z); iden := -sin(-z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("sin(z) = sin(z+2Pi) periodic"); for i in p'range loop begin z := p(i); base := sin(z); iden := sin(z+two_pi); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("conjugate(sin(z)) = sin(conjugate(z))"); for i in p'range loop begin z := p(i); base := conjugate(sin(z)); iden := sin(conjugate(z)); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("test cos"); put_line("cos(arccos(z))=z"); for i in p'range loop begin z := p(i); base := z; iden := cos(arccos(z)); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("cos(z)=(exp(iz)+exp(-iz))/2"); for i in p'range loop begin z := p(i); base := cos(z); iden := (exp(complex_i*z)+exp(-complex_i*z))/2.0; track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("cos(z) = cosh iz"); for i in p'range loop begin z := p(i); base := cos(z); iden := cosh(complex_i*z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("cos(z) = sin(z + Pi/2)"); for i in p'range loop begin z := p(i); base := cos(z); iden := sin(z + Pi_2); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("cos(z) = cos(-z)"); for i in p'range loop begin z := p(i); base := cos(z); iden := cos(-z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("cos(z) = cos(z+2 Pi) periodic"); for i in p'range loop begin z := p(i); base := cos(z); iden := cos(z+two_pi); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("conjugate(cos(z)) = cos(conjugate(z))"); for i in p'range loop begin z := p(i); base := conjugate(cos(z)); iden := cos(conjugate(z)); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("test tan"); put_line("tan(arctan(z))=z"); for i in p'range loop begin z := p(i); base := z; iden := tan(arctan(z)); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("tan(z)=-i(exp(iz)-exp(-iz))/(exp(iz)+exp(-iz))"); for i in p'range loop begin z := p(i); base := tan(z); iden := -complex_i*(exp(complex_i*z)-exp(-complex_i*z))/ (exp(complex_i*z)+exp(-complex_i*z)); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("tan(z) = sin x cos x/(cos2 x+ sinh2 x) +i sinh y cosh y/(...)"); for i in p'range loop declare x, y , denom : real; begin z := p(i); x := RE(z); y := IM(z); base := tan(z); begin denom := rf.cos(x)*rf.cos(x)+rf.sinh(y)*rf.sinh(y); iden := compose_from_cartesian(rf.sin(x)*rf.cos(x)/denom, rf.sinh(y)*rf.cosh(y)/denom); exception when constraint_error | numeric_error => iden := compose_from_cartesian(0.0,1.0); if y<0.0 then iden := compose_from_cartesian(0.0,-1.0); end if; end; track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("tan(z) = -i tanh iz"); for i in p'range loop begin z := p(i); base := tan(z); iden := -complex_i*tanh(complex_i*z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("tan(z) = cot(Pi/2-z)"); for i in p'range loop begin z := p(i); base := tan(z); iden := cot(pi_2-z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("tan(z) = 1/cot(z)"); for i in p'range loop begin z := p(i); base := tan(z); iden := 1.0/cot(z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("tan(z) = -tan(-z)"); for i in p'range loop begin z := p(i); base := tan(z); iden := -tan(-z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("tan(z) = tan(z+2 Pi) periodic"); for i in p'range loop begin z := p(i); base := tan(z); iden := -tan(-z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("conjugate(tan(z)) = tan(conjugate(z))"); for i in p'range loop begin z := p(i); base := conjugate(tan(z)); iden := tan(conjugate(z)); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("test cot"); put_line("cot(z)=i(exp(iz)+exp(-iz))/(exp(iz)-exp(-iz))"); for i in p'range loop begin z := p(i); base := cot(z); iden := complex_i*(exp(complex_i*z)+exp(-complex_i*z))/ (exp(complex_i*z)-exp(-complex_i*z)); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("cot(z) = sin x cos x/(sin2 x+ sinh2 x) -i sinh y cosh y/(...)"); for i in p'range loop declare x, y , denom : real; begin z := p(i); x := RE(z); y := IM(z); base := cot(z); begin denom := rf.sin(x)*rf.sin(x)+rf.sinh(y)*rf.sinh(y); iden := compose_from_cartesian(rf.sin(x)*rf.cos(x)/denom, -rf.sinh(y)*rf.cosh(y)/denom); exception when constraint_error | numeric_error => iden := compose_from_cartesian(0.0,-1.0); if y<0.0 then iden := compose_from_cartesian(0.0,1.0); end if; end; track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("cot(z)=1/tan(z)"); for i in p'range loop begin z := p(i); base := cot(z); iden := 1.0/tan(z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("cot(z) = tan(pi/2-z)"); for i in p'range loop begin z := p(i); base := cot(z); iden := tan(pi_2-z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("cot(z) = -cot(-z)"); for i in p'range loop begin z := p(i); base := cot(z); iden := -cot(-z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("cot(z) = cot(z+2 Pi) periodic"); for i in p'range loop begin z := p(i); base := cot(z); iden := cot(z+two_pi); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("conjugate(cot(z)) = cot(conjugate(z))"); for i in p'range loop begin z := p(i); base := conjugate(cot(z)); iden := cot(conjugate(z)); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("test arcsin"); put_line("arcsin(sin(z)) = z real part in [-Pi/2,Pi/2]"); for i in p'range loop begin z := p(i); if in_defined_range_arcsin(z) and then in_principal_domain_sin(z) and then not on_branch_cut_arcsin(sin(z)) then base := z; iden := arcsin(sin(z)); track(z,base,iden); else excl := excl + 1; end if; exception when others => excp := excp + 1; end; end loop; print_track; put_line("arcsin(z)) = pi/2 - arccos(z)"); for i in p'range loop begin z := p(i); base := arcsin(z); iden := pi_2-arccos(z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("arcsin(z)) = -i arcsinh(iz)"); for i in p'range loop begin z := p(i); base := arcsin(z); iden := -complex_i*arcsinh(complex_i*z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("arcsin(z) = -i log(iz+sqrt(1-z)*sqrt(1+z))"); for i in p'range loop begin z := p(i); base := arcsin(z); iden := -complex_i*log(complex_i*z+sqrt(complex_one-z)* sqrt(complex_one+z)); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("arcsin(z) = real part +i imag part"); for i in p'range loop begin z := p(i); base := arcsin(z); iden := local_arcsin(z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("arcsin(z)) = arctan(z/sqrt(1-z*z)) not on any branch cuts"); for i in p'range loop begin if not on_branch_cut_arcsin(z) and not on_branch_cut_arctan(z) and not on_branch_cut_sqrt(z) then z := p(i); base := arcsin(z); iden := arctan(z/sqrt((1.0-z)*(1.0+z))); track(z,base,iden); else excl := excl + 1; end if; exception when others => excp := excp + 1; end; end loop; print_track; put_line("conjugate(arcsin(z)) = arcsin(conjugate(z)) not on branch cut"); for i in p'range loop begin z := p(i); if not on_branch_cut_arcsin(z) then base := conjugate(arcsin(z)); iden := arcsin(conjugate(z)); track(z,base,iden); else excl := excl + 1; end if; exception when others => excp := excp + 1; end; end loop; print_track; put_line("test arccos"); put_line("arccos(z) = -2i ln(sqrt((1+z)/2) +i sqrt((1-z)/2))"); for i in p'range loop begin z := p(i); base := arccos(z); iden := -complex_i*2.0*log(sqrt((1.0+z)/2.0) + complex_i*sqrt((1.0-z)/2.0)); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("arccos(z) = -i ln(z+i sqrt(1+z)*sqrt(1-z)))"); for i in p'range loop begin z := p(i); base := arccos(z); iden := -complex_i*log(z + complex_i*sqrt(1.0+z)* sqrt(1.0-z)); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("arccos(z) = Pi/2 - arcsin(z)"); for i in p'range loop begin z := p(i); base := arccos(z); iden := pi_2-arcsin(z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("arccos(z) = real part +i imag part"); for i in p'range loop begin z := p(i); base := arccos(z); iden := local_arccos(z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("arccos(cos(z)) = z not on branch cut, principal domain"); for i in p'range loop begin z := p(i); if in_defined_range_arccos(z) and then in_principal_domain_cos(z) and then not on_branch_cut_arccos(cos(z)) then base := z; iden := arccos(cos(z)); track(z,base,iden); end if; exception when others => excp := excp + 1; end; end loop; print_track; put_line("arccos(z)) = arctan(sqrt(1-z*z)/z) not on any branch cuts, adjust"); for i in p'range loop begin if not on_branch_cut_arccos(z) and not on_branch_cut_arctan(z) and not on_branch_cut_sqrt(z) then z := p(i); base := arccos(z); iden := arctan(sqrt((1.0-z)*(1.0+z))/z); track(z,base,iden); else excl := excl + 1; end if; exception when others => excp := excp + 1; end; end loop; print_track; put_line("conjugate(arccos(z)) = arccos(conjugate(z))"); for i in p'range loop begin z := p(i); if not on_branch_cut_arccos(z) then base := conjugate(arccos(z)); iden := arccos(conjugate(z)); track(z,base,iden); else excl := excl + 1; end if; exception when others => excp := excp + 1; end; end loop; print_track; put_line("test arctan"); put_line("arctan(sin(z)/cos(z)) = z"); for i in p'range loop begin z := p(i); if RE(z) > -pi_2 and RE(z) < pi_2 then base := z; iden := arctan(sin(z)/cos(z)); track(z,base,iden); end if; exception when others => excp := excp + 1; end; end loop; print_track; put_line("arctan(z) = pi/2 -arccot(1/z)"); for i in p'range loop begin z := p(i); base := arctan(z); iden := pi_2-arccot(z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("arctan(z) = arccot(1/z) not on either branch cut, adjust domain"); for i in p'range loop begin z := p(i); if not on_branch_cut_arctan(z) and not on_branch_cut_arccot(1.0/z) then base := arctan(z); iden := arccot(1.0/z); if RE(iden) > pi_2 then iden := iden - pi; end if; track(z,base,iden); else excl := excl + 1; end if; exception when others => excp := excp + 1; end; end loop; print_track; put_line("arctan(z) = -i arctanh iz"); for i in p'range loop begin z := p(i); base := arctan(z); iden := -complex_i*arctanh(complex_i*z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("conjugate(arctan(z)) = arctan(conjugate(z)) not on branch cut"); for i in p'range loop begin z := p(i); if not on_branch_cut_arctan(z) then base := conjugate(arctan(z)); iden := arctan(conjugate(z)); track(z,base,iden); else excl := excl + 1; end if; exception when others => excp := excp + 1; end; end loop; print_track; put_line("test arccot"); put_line("arccot(cos(z)/sin(z)) = z"); for i in p'range loop begin z := p(i); if RE(z) > -pi_2 and RE(z) < pi_2 then base := z; iden := arccot(cos(z)/sin(z)); track(z,base,iden); end if; exception when others => excp := excp + 1; end; end loop; print_track; put_line("arccot z = Pi/2 - arctan z"); for i in p'range loop begin z := p(i); base := arccot(z); iden := pi_2-arctan(z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("arccot z = arctan 1/z"); for i in p'range loop begin z := p(i); base := arccot(z); iden := arctan(1.0/z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("conjugate(arccot(z)) = arccot(conjugate(z))"); for i in p'range loop begin z := p(i); base := conjugate(arccot(z)); iden := arccot(conjugate(z)); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("test sinh"); put_line("sinh(arcsinh(z))=z"); for i in p'range loop begin z := p(i); base := z; iden := sinh(arcsinh(z)); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("sinh(z)=(exp(z)-exp(-z))/2"); for i in p'range loop begin z := p(i); base := sinh(z); iden := (exp(z)-exp(-z))/(0.0,2.0); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("sinh(z)=(exp(z)-exp(-z))/2"); for i in p'range loop begin z := p(i); base := sinh(z); iden := (exp(z)-exp(-z))/2.0; track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("sinh(z)=-i sin(iz)"); for i in p'range loop begin z := p(i); base := sinh(z); iden := -complex_i*sin(complex_i*z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("sinh(z) = -i cosh(z +i Pi/2)"); for i in p'range loop begin z := p(i); base := sinh(z); iden := -complex_i*cosh(z+complex_i*Pi_2); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("sinh(z) = -sinh(-z) odd function"); for i in p'range loop begin z := p(i); base := sinh(z); iden := -sinh(-z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("sinh(z)=sinh(z + i 2Pi) periodic"); for i in p'range loop begin z := p(i); base := sinh(z); iden := sinh(z+complex_i*two_pi); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("conjugate(sinh(z)) = sinh(conjugate(z))"); for i in p'range loop begin z := p(i); base := conjugate(sinh(z)); iden := sinh(conjugate(z)); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("test cosh"); put_line("cosh(arccosh(z))=z"); for i in p'range loop begin z := p(i); base := z; iden := cosh(arccosh(z)); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("cosh(z)=(exp(z)+exp(-z))/2.0"); for i in p'range loop begin z := p(i); base := cosh(z); iden := (exp(z)+exp(-z))/2.0; track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("cosh(z)=cos(i z)"); for i in p'range loop begin z := p(i); base := cosh(z); iden := cos(complex_i*z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("cosh(z) = -i sinh(z +i Pi/2)"); for i in p'range loop begin z := p(i); base := cosh(z); iden := -complex_i*sinh(z+complex_i*Pi_2); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("cosh(z) = cosh(-z)"); for i in p'range loop begin z := p(i); base := cosh(z); iden := cosh(-z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("cosh(z)=cosh(z + i 2Pi) periodic"); for i in p'range loop begin z := p(i); base := cosh(z); iden := cosh(z+complex_i*two_pi); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("conjugate(cosh(z)) = cosh(conjugate(z))"); for i in p'range loop begin z := p(i); base := conjugate(cosh(z)); iden := cosh(conjugate(z)); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("test tanh"); put_line("tanh(arctanh(z))=z"); for i in p'range loop begin z := p(i); base := z; iden := tanh(arctanh(z)); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("tanh(z)=(exp(z)-exp(-z))/(exp(z)+exp(-z))"); for i in p'range loop begin z := p(i); base := tanh(z); iden := (exp(z)-exp(-z))/(exp(z)+exp(-z)); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("tanh(z) = sinh x cosh x/(sinh2 x+ cos2 x) +i sin y cos y/(...)"); for i in p'range loop declare x, y , denom : real; begin z := p(i); x := RE(z); y := IM(z); base := tanh(z); begin denom := rf.sinh(x)*rf.sinh(x)+rf.cos(y)*rf.cos(y); iden := compose_from_cartesian(rf.sinh(x)*rf.cosh(x)/denom, rf.sin(y)*rf.cos(y)/denom); exception when constraint_error | numeric_error => iden := compose_from_cartesian(0.0,1.0); if y<0.0 then iden := compose_from_cartesian(0.0,-1.0); end if; end; track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("tanh(z) = -i tan iz"); for i in p'range loop begin z := p(i); base := tanh(z); iden := -complex_i*tan(complex_i*z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("tanh(z) = coth(z +i Pi/2"); for i in p'range loop begin z := p(i); base := tanh(z); iden := coth(z+complex_i*Pi_2); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("conjugate(tanh(z)) = tanh(conjugate(z))"); for i in p'range loop begin z := p(i); base := conjugate(tanh(z)); iden := tanh(conjugate(z)); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("test coth"); put_line("coth(z) = tanh(z +i Pi/2)"); for i in p'range loop begin z := p(i); base := coth(z); iden := tanh(z+complex_i*Pi_2); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("coth(z)=(exp(z)+exp(-z))/(exp(z)-exp(-z))"); for i in p'range loop begin z := p(i); base := coth(z); iden := (exp(z)+exp(-z))/(exp(z)-exp(-z)); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("coth(z) = sinh x cosh x/(sinh2 x+ sin2 y) +i sin y cos y/(...)"); for i in p'range loop declare x, y , denom : real; begin z := p(i); x := RE(z); y := IM(z); base := coth(z); begin denom := rf.sinh(x)*rf.sinh(x)+rf.sin(y)*rf.sin(y); iden := compose_from_cartesian(rf.sinh(x)*rf.cosh(x)/denom, -rf.sin(y)*rf.cos(y)/denom); -- exception -- when constraint_error | numeric_error => -- iden := compose_from_cartesian(0.0,1.0); -- if y<0.0 then -- iden := compose_from_cartesian(0.0,-1.0); wrong for (-1000,-1000) -- end if; end; track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("coth(z) = i cot(iz))"); for i in p'range loop begin z := p(i); base := coth(z); iden := complex_i*cot(complex_i*z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("conjugate(coth(z)) = coth(conjugate(z))"); for i in p'range loop begin z := p(i); base := conjugate(coth(z)); iden := coth(conjugate(z)); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("test arcsinh"); put_line("conjugate(arcsinh(z)) = arcsinh(conjugate(z))"); for i in p'range loop begin z := p(i); if IM(z) > -pi_2 and IM(z) < pi_2 then base := conjugate(arcsinh(z)); iden := arcsinh(conjugate(z)); track(z,base,iden); end if; exception when others => excp := excp + 1; end; end loop; print_track; put_line("test arccosh"); put_line("conjugate(arccosh(z)) = arccosh(conjugate(z))"); for i in p'range loop begin z := p(i); if RE(z) > 0.0 and IM(z) > -pi_2 and IM(z) < pi_2 then base := conjugate(arccosh(z)); iden := arccosh(conjugate(z)); track(z,base,iden); end if; exception when others => excp := excp + 1; end; end loop; print_track; put_line("test arctanh"); put_line("arctanh(tanh(z)) = z not on branch cut, principal domain"); for i in p'range loop begin z := p(i); if in_principal_domain_tanh(z) and then not on_branch_cut_arctanh(tanh(z)) then base := z; iden := arctanh(tanh(z)); track(z,base,iden); else excl := excl + 1; end if; exception when others => excp := excp + 1; end; end loop; print_track; put_line("arctanh(z) = -i arctan(iz)"); for i in p'range loop begin z := p(i); base := arctanh(z); iden := -complex_i*arctan(complex_i*z); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("arctanh(sinh(z)/cosh(z)) = z not on branch cut, principal"); for i in p'range loop begin z := p(i); if in_principal_domain_sinh(z) and then not on_branch_cut_arctanh(sinh(z)/cosh(z)) then base := z; iden := arctanh(sinh(z)/cosh(z)); track(z,base,iden); else excl := excl + 1; end if; exception when others => excp := excp + 1; end; end loop; print_track; put_line("arctanh(z) = arccoth(z) - i Pi/2 not on either branch cut"); for i in p'range loop begin z := p(i); if not on_branch_cut_arctanh(z) and not on_branch_cut_arccoth(z) then base := arctanh(z); iden := arccoth(z); track(z,base,iden); else excl := excl + 1; end if; exception when others => excp := excp + 1; end; end loop; print_track; put_line("conjugate(arctanh(z)) = arctanh(conjugate(z)) not on branch cut"); for i in p'range loop begin z := p(i); if not on_branch_cut_arctanh(z) then base := conjugate(arctanh(z)); iden := arctanh(conjugate(z)); track(z,base,iden); else excl := excl + 1; end if; exception when others => excp := excp + 1; end; end loop; print_track; put_line("arccoth"); put_line("arccoth(coth(z)) = z not on branch cut, principal domain"); for i in p'range loop begin z := p(i); if not on_branch_cut_arccoth(z) and in_principal_domain_coth(z) then base := z; iden := arccoth(coth(z)); track(z,base,iden); else excl := excl + 1; end if; exception when others => excp := excp + 1; end; end loop; print_track; put_line("arccoth z = 1/2(log(z+1)-log(z-1))"); for i in p'range loop begin z := p(i); base := arccoth(z); iden := 0.5*(log(z+1.0)-log(z-1.0)); track(z,base,iden); exception when others => excp := excp + 1; end; end loop; print_track; put_line("arccoth z = arctanh z +i Pi/2 not on either branch cut"); for i in p'range loop begin z := p(i); if not on_branch_cut_arccoth(z) and not on_branch_cut_arctanh(z) then base := arccoth(z); iden := arctanh(z)+complex_i*Pi_2; track(z,base,iden); else excl := excl + 1; end if; exception when others => excp := excp + 1; end; end loop; print_track; put_line("arccoth z = arctanh 1/z offset range, not on either branch cut"); for i in p'range loop begin z := p(i); if not on_branch_cut_arccoth(z) and then not on_branch_cut_arctanh(1.0/z) then base := arccoth(z); iden := arctanh(1.0/z); if IM(iden) < 0.0 then -- adjust range iden := iden + complex_i*pi; end if; track(z,base,iden); else excl := excl + 1; end if; exception when others => excp := excp + 1; end; end loop; print_track; put_line("conjugate(arccoth(z)) = arccoth(conjugate(z)) not on branch cut," & " adjust range"); for i in p'range loop begin z := p(i); if not on_branch_cut_arccoth(z) then base := conjugate(arccoth(z)); if IM(base) < 0.0 then base := base + complex_i*pi; end if; iden := arccoth(conjugate(z)); track(z,base,iden); else excl := excl + 1; end if; exception when others => excp := excp + 1; end; end loop; print_track; put_line("end test complex identities"); exception when others => put_line("abnormal exception, end test complex identities"); end test_generic_complex_identities; -- instantiations to produce procedures that can be linked and executed with test_generic_complex_identities; procedure test_complex_identities is new test_generic_complex_identities(float); with test_generic_complex_identities; procedure test_long_complex_identities is new test_generic_complex_identities(long_float);