!$Author$ !$Date$ !$Revision$ !$HeadURL$ double precision function enorm(n,x) integer n double precision x(n) ! ********** ! ! function enorm ! ! given an n-vector x, this function calculates the ! euclidean norm of x. ! ! the euclidean norm is computed by accumulating the sum of ! squares in three different sums. the sums of squares for the ! small and large components are scaled so that no overflows ! occur. non-destructive underflows are permitted. underflows ! and overflows do not occur in the computation of the unscaled ! sum of squares for the intermediate components. ! the definitions of small, intermediate and large components ! depend on two constants, rdwarf and rgiant. the main ! restrictions on these constants are that rdwarf**2 not ! underflow and rgiant**2 not overflow. the constants ! given here are suitable for every known computer. ! ! the function statement is ! ! double precision function enorm(n,x) ! ! where ! ! n is a positive integer input variable. ! ! x is an input array of length n. ! ! subprograms called ! ! fortran-supplied ... dabs,dsqrt ! ! argonne national laboratory. minpack project. march 1980. ! burton s. garbow, kenneth e. hillstrom, jorge j. more ! ! ********** integer i double precision agiant,floatn,one,rdwarf,rgiant,s1,s2,s3,xabs, & & x1max,x3max,zero data one,zero,rdwarf,rgiant /1.0d0,0.0d0,3.834d-20,1.304d19/ s1 = zero s2 = zero s3 = zero x1max = zero x3max = zero floatn = n agiant = rgiant/floatn do 90 i = 1, n xabs = dabs(x(i)) if (xabs .gt. rdwarf .and. xabs .lt. agiant) go to 70 if (xabs .le. rdwarf) go to 30 ! ! sum for large components. ! if (xabs .le. x1max) go to 10 s1 = one + s1*(x1max/xabs)**2 x1max = xabs go to 20 10 continue s1 = s1 + (xabs/x1max)**2 20 continue go to 60 30 continue ! ! sum for small components. ! if (xabs .le. x3max) go to 40 s3 = one + s3*(x3max/xabs)**2 x3max = xabs go to 50 40 continue if (xabs .ne. zero) s3 = s3 + (xabs/x3max)**2 50 continue 60 continue go to 80 70 continue ! ! sum for intermediate components. ! s2 = s2 + xabs**2 80 continue 90 continue ! ! calculation of norm. ! if (s1 .eq. zero) go to 100 enorm = x1max*dsqrt(s1+(s2/x1max)/x1max) go to 130 100 continue if (s2 .eq. zero) go to 110 if (s2 .ge. x3max) & & enorm = dsqrt(s2*(one+(x3max/s2)*(x3max*s3))) if (s2 .lt. x3max) & & enorm = dsqrt(x3max*((s2/x3max)+(x3max*s3))) go to 120 110 continue enorm = x3max*dsqrt(s3) 120 continue 130 continue return ! ! last card of function enorm. ! end