#C19.txt: ExpMath(RU)
Help:=proc(): print(` SofS(n) , P(L), Q(L), FibS(n) , wt(s,a) `): end:

#SofS(n): the set of all the pairs [a,b], a<=b,  such that a^2+b^2=n
SofS:=proc(n) local S,a,b:

S:={}:

for a from 1 to trunc(sqrt(n))+1 do
  b:=sqrt(n-a^2):
  if type(b,integer)  and a<=b then
    S:=S union {[a,b]}:
  fi:
od:

S:

end:



#The numerator of the (simple) continued fraction L
P:=proc(L) option remember:

if L=[] then
 1:
elif nops(L)=1 then
  L[1]:
else
expand(L[1]*P([op(2..nops(L),L)])+ P([op(3..nops(L),L)])):

fi:
 end:


Q:=proc(L):
P([op(2..nops(L),L)]):
end:

#FibS(n): the set of vectors in {1,2} that add-up to n
FibS:=proc(n) local S1,S2,s:
option remember:
if n=0 then
 RETURN({[]}):
elif n=1 then
  RETURN({[1]}):
else
 S1:=FibS(n-1):
 S2:=FibS(n-2):
  
 {seq([1,op(s)],s in S1), seq([2,op(s)],s in S2)}:
fi:
end:

#wt(s): the weight of the Fibonacci word s, product of a[j] over all intermediate stops
#j that follow immediately j-1
#wt(122212)=a[1]*a[8] wt(22222)=1
wt:=proc(s,a) local i,p:
p:=1:
for i from 1 to nops(s) do
 if s[i]=1 then
  p:=p*a[add(s[j],j=1..i)]:
 fi:
od:

p:

end:

#CheckQI(L,k): Checks the identity P([a1,...., an])=P([a1...,ak])*P([a(k+1)... an])+  ... 
CheckQI:=proc(L,k) local i:
expand(P(L)-P([op(1..k,L)])*P([op(k+1..nops(L),L)])-
P([op(1..k-1,L)])*P([op(k+2..nops(L),L)])):
end:






###from C18.txt
#C18.txt March 26, 2018, ExpMath (RU)
Help18:=proc(): print(`RatToCF(x),NuToCF(x,n) , EvalCF(L), SymCF(a,n) `):end:

#RatToCF(x): converts a pos. rational number x into a continued fraction
RatToCF:=proc(x) local a,x1:
 if not (x>0 and type(x,rational)) then
  RETURN(FAIL):
 fi:

a:=trunc(x):

if a=x then
  RETURN([a]):
fi:

x1:=1/(x-a):

[a, op(RatToCF(x1))]:


end:


#NuToCF(x,n): appx. a pos. number x into a continued fraction with <=n terms
NuToCF:=proc(x,n) local a,x1:

 
 if not (x>0 and type(n,integer) and n>=1) then
  RETURN(FAIL):
 fi:


a:=trunc(x):

if n=1 then
  RETURN([a]):
fi:
  
if x=a then
  RETURN([a]):
fi:

x1:=1/(x-a):

[a,op(NuToCF(x1,n-1))]:
end:

#EvalCF(L): evaluates the simple continued 
EvalCF:=proc(L) local L1,x:
if not (type(L,list) and nops(L)>0) then
  RETURN(FAIL):
fi:

if nops(L)=1 then
 RETURN(L[1]):
fi:

L1:=[op(2..nops(L),L)]:

x:=EvalCF(L1):

normal(L[1]+1/x):
end:

#SymCF: the rational function of a[1], ..., a[n] describing EvalCF([a[1], ..., a[n]]):
#the numerator followed by the denominator
SymCF:=proc(a,n) local i,yonah:
yonah:=EvalCF([seq(a[i],i=1..n)]):

numer(yonah),denom(yonah):
end:



####