subroutine distancesys(nunknowns,scldxp,
     &  scldfequ,scldfsqsum,idowhat)
      implicit none
!
!--------------------- Variables through arguments -------------------------
!(in) nunknowns:          the number of unknowns and non-linear equations
!(in) scldxp(1:nunknowns):  the scaled unknowns in the nonlinear system
!(out) scldfequ(1:nunknowns): the scaled function values evaluated at imported unknowns
!(out) scldfsqsum: half of the sum of the squared scaled function values.
      integer nunknowns,idowhat
      double precision scldxp(1:nunknowns),
     &  scldfequ(1:nunknowns),scldfsqsum
!--------------------------- Local variables-----------------------------------
      integer i,j,numys
      parameter
      double precision xvar(nunknowns),xscalingfact(nunknowns),
     &  xtarget(nunknowns),yvar(numys),yscalingfact(numys),
     &  ytarget(numys),diffy(numys)
!^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
      do i=1,nunknowns
        scalingfact(i)=1.0d0
        xvar(i)=scldxp(i)/xscalingfact(i)
      enddo
      call surface(numys,yvar,nunknowns,xvar)
      do i=1,numys
        diff(i)=yvar(i)-ytarget(i)
      endif
      if(idowhat.eq.1)then
        scldfsqsum=0.0d0
        do i=1,nunknowns
          scldfequ(i)=scldfequ(i)*scalingfact(i)
          scldfsqsum=scldfsqsum+scldfequ(i)*scldfequ(i)        
        enddo  
        scldfsqsum=0.5d0*scldfsqsum
      endif
      if(idowhat.eq.2)then
        
      return
      end subroutine leafsys
      
      subroutine surface(numys,ysurface,nunknowns,xvar)
      implicit none
      integer numys,nunknowns
      double precision ysurface(numys),xvar(nunknowns)
      ysurface(1)=xvar(1)
      return
      end

!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
      subroutine sqsum_fequ(nunknowns,xp,fequsqsum)
!This subroutine returns half of the sum of the squared equation residues 
      implicit none
      integer nunknowns
      double precision xp(nunknowns),fequsqsum,
     &  fequ(nunknowns)
      integer idowhat
      parameter(idowhat=1)
      call distancesys(nunknowns,xp,fequ,fequsqsum,idowhat)
      return
      end

      double precision function cpf1dim_fequ(x)
!this function subroutine returns half of the sum of the squared equation residues for line search
      INTEGER NMAX
      double precision x
      PARAMETER (NMAX=1000)
      INTEGER j,ncom,idowhat
      parameter(idowhat=1)
      double precision pcom(NMAX),xicom(NMAX),
     &  xt(NMAX),fequ(NMAX)
      COMMON /cpf1com/ pcom,xicom,ncom
      save /cpf1com/
      do 11 j=1,ncom
        xt(j)=pcom(j)+x*xicom(j)
11    continue
      call distancesys(ncom,xt,fequ,cpf1dim_fequ,idowhat)
      return
      END
      
      subroutine getequationvalues((nunknowns,xp,
     &  fequ,fequsqsum)
!this subroutine is for solving a nonlinear system
      implicit none
      integer nunknowns,idowhat
      double precision xp(nunknowns),fequ(nunknowns),fequsqsum
      parameter(idowhat=1)
      call distancesys(nunknowns,xp,fequ,fequsqsum,idowhat)
      return
      end

      subroutine getsqdistance(nunknowns,xp,sqdistance)
!This subroutine returns the distance between a point on a curve (surface) and another point 
!that is specified in distancesys 
      implicit none
      integer nunknowns,idowhat
      double precision xp(nunknowns),sqdistance,fequ(nunknowns)
      parameter(idowhat=2)
      call distancesys(nunknowns,xp,fequ,sqdistance,idowhat)
      return
      end

      double precision function cpf1dim_distance(x)
!this function subroutine returns the distance for line search
      INTEGER NMAX
      double precision x
      PARAMETER (NMAX=1000)
      INTEGER j,ncom,idowhat
      parameter(idowhat=2)
      double precision pcom(NMAX),xicom(NMAX),
     &  xt(NMAX),fequ(NMAX)
      COMMON /cpf1com/ pcom,xicom,ncom
      save /cpf1com/
      do 11 j=1,ncom
        xt(j)=pcom(j)+x*xicom(j)
11    continue
      call distancesys(ncom,xt,fequ,cpf1dim_distance,idowhat)
      return
      END
!@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@