> 
f:=array([seq(evalf(2*x[i]+x[i]^2+0.99/i^2),i=1..N)]);

> 
fsolve({f[1],f[2],f[3],f[4]},{x[1]=1.2,x[2]=1,x[3]=0,x[4]=0});

> 
eqsA := [seq(F[i]=f[i],i=1..N)]:

> 
irform2 := StatSeq(seq(Assign(op(i)),i=eqsA)):

> 
prccons:= codegen[intrep2maple](Proc(Parameters(x::Array,F::Array),irform2)):

> 
f3:=Compiler:Compile(prccons):

> 
eqsJ := [seq(seq(Jac[i,j]=evalf(diff(f[i],x[j])),j=1..N),i=1..N)];

> 
irformJ := StatSeq(seq(Assign(op(i)),i=eqsJ)):

> 
prcconsJ:= codegen[intrep2maple](Proc(Parameters(x::Array,Jac::Array),irformJ)):

> 
J3:=Compiler:Compile(prcconsJ):

Following linear sovler algorithm is basiclly from dsolve/numeric/SC
.
> 
s3:=proc(n::posint, A::Array( datatype = float[ 8 ] ) , ip::Array( datatype = integer[ 4 ] ),b::Array( datatype = float[ 8 ] ) ) local i::integer, j::integer, k::integer, m::integer, t::float; ip[n] := 1; for k to n1 do m := k; for i from k+1 to n do if abs(A[m,k]) < abs(A[i,k]) then m := i end if end do; ip[k] := m; if m <> k then ip[n] := ip[n] end if; t := A[m,k]; A[m,k] := A[k,k]; A[k,k] := t; if t = 0 then ip[n] := 0; #return ip end if; for i from k+1 to n do A[i,k] := A[i,k]/t end do; for j from k+1 to n do t := A[m,j]; A[m,j] := A[k,j]; A[k,j] := t; if t <> 0 then for i from k+1 to n do A[i,j] := A[i,j]+A[i,k]*t end do end if end do end do; if A[n,n] = 0 then ip[n] := 0 end if; #ip[n] if ip[n] = 0 then error "The matrix A is numerically singular" end if; if 1 < n then for k to n1 do m := ip[k]; t := b[m]; b[m] := b[k]; b[k] := t; for i from k+1 to n do b[i] := b[i]+A[i,k]*t end do end do; for k from n by 1 to 2 do b[k] := b[k]/A[k,k]; t := b[k]; for i to k1 do b[i] := b[i]+A[i,k]*t end do end do end if; b[1] := b[1]/A[1,1]; return 0.0; end proc:

> 
s3c:=Compiler:Compile(s3):

> 
Newton:=proc(x::Array(datatype=float[8]), f0::Array(datatype=float[8]), j0::Array(datatype=float[8]), N::integer, ip::Array(datatype=integer[4])) global f3,J3,s3c; local err::float[8],i::integer,k::integer,c::integer; err:=10.0; c:=0; while err>1e8 and c<100 do f3(x,f0); J3(x,j0); for i from 1 to N do for k from 1 to N do j0[i,k]:=1.0*j0[i,k]; end do; end do; s3c(N,j0,ip,f0); for i from 1 to N do x[i] := x[i]+f0[i]; end do; f0[1] := abs(f0[1]); err := abs(f0[1]); for i from 2 to N do f0[i] := abs(f0[i]); if f0[i]>err then err := f0[i]; end if; end do; err:=err/N; c:=c+1; end do: err; end proc:

> 
Newtonc := Compiler:Compile(Newton):

Warning, the function names {J3, f3, s3c} are not recognized in the target language
> 
x0:=Array(1..N,datatype=float[8]): f0:=Array(1..N,datatype=float[8]): j0:=Array(1..N,1..N,datatype=float[8]): ip:=Array(1..N,datatype=integer[4]):

> 
ArrayTools:Fill(0.5,x0): ArrayTools:Fill(0.0,f0): ArrayTools:Fill(0.0,j0): ArrayTools:Fill(1,ip): x0; f0; j0, ip; t11:=time(): for ii from 1 to 10000 do Newtonc(x0,f0,j0,N,ip); # compiled end do: time()t11; x0, f0, j0, ip;

> 
ArrayTools:Fill(0.5,x0): ArrayTools:Fill(0.0,f0): ArrayTools:Fill(0.0,j0): ArrayTools:Fill(1,ip): x0; f0; j0, ip; t11:=time(): for ii from 1 to 10000 do Newton(x0,f0,j0,N,ip); #uncompiled end do: time()t11; x0, f0, j0, ip;

