# ClosedForm.ml

```(* Copyright INRIA and Microsoft Corporation, 2008-2013. *)

INCLUDE "preamble.ml"

let title _ = <:text<Closed Form>>

(* A wrapper function to call DDMF:-Hypergeometric:-generalTerms *)
(* that does some pre- and post-processing. *)
let call_generalTerms recsys inits u n x lps =
(* Write recsys in terms of positive shifts (required by generalTerms). *)
let shift =
<< min(map(a -> op(1, a) - \$(n), indets(\$(recsys), function))) >> in
let recsys = << subs(\$(n) = \$(n) - \$(shift), \$(recsys)) >> in
(* Compute closed forms of the unknown sequences u[i](n). *)
let p = << indets([\$(recsys), \$(inits)], name) minus {\$(n)} >> in
let general_terms =
<<
DDMF:-Hypergeometric:-generalTerms(\$(recsys), \$(inits), \$(n), \$(x), \$(p))
>> in
<< [seq(\$(general_terms)[\$(u)[i]], i in \$(lps))] >>

(* Compute a closed form expression of a transseries record. The closed form *)
(* is returned as a Maple expression. The record itself is updated, when the *)
(* closed form is computed for the first time. *)
let closed_form_of_transseries ser =
(* Shorter name for the constituents of the transseries record. *)
let x = << \$(ser):-variable >>
and ram = << \$(ser):-ramification >>
and u = << \$(ser):-sequence_name >>
and n = << \$(ser):-index_name >> in

(* A local variable to ensure that the final expression is displayed in *)
(* terms of ser:-variable, i.e., the "transf" given in the sf database, *)
(* and to make intermediate operations simpler. *)
let x_loc = << proc() local x_loc; x_loc end proc() >> in

(* Whether the power factor x^alpha is displayed inside the sum or not. *)
let pow_factor, x_expon =
if <:bool< type(\$(ram) * \$(ser):-alpha, integer) >>
then << 1 >>, << \$(ser):-alpha >>
else << combine(\$(x) ^ \$(ser):-alpha, 'power', 'symbolic') >>, << 0 >> in

(* Maple expression of an element in the list :-general. *)
let expr_of_general_term gen =
let x_part = << \$(x_loc) ^ (\$(gen):-period / \$(ram) + \$(x_expon)) >>
and cf = << \$(gen):-cform >> in
let facs = << `if`(op(0, \$(cf)) = `*`, [op(\$(cf))], [\$(cf)]) >> in
let facs = << [selectremove(has, \$(facs), \$(n))] >> in
let out_fac = << `*`(op(\$(facs)[2])) >>
and in_fac =
<< map(a -> sort(a, indets(a, name), descending), \$(facs)[1]) >> in
let in_fac = << __DynaMoW_times([op(sort(\$(in_fac))), \$(x_part)]) >> in
let in_fac = << DDMF:-dynamow_times_to_frac(\$(in_fac), []) >> in
<< \$(out_fac) * Sum(\$(in_fac), \$(n) = 0 .. infinity) >>
in

(* Update the transseries record, concerning the closed form sol which *)
(* is related to index i (which corresponds to u[i](n) * log(x)^i). *)
let update_index i sol =
let ln = << ln(\$(x_loc))^\$(i) >> in
if <:bool< \$(sol):-isHyper >> then
(* If a closed form is found, it is stored in the transseries record. *)
let gen = CommonTools.symb_list_of_symb << \$(sol):-general >> in
let exprs = List.map expr_of_general_term gen in
let expr = << `+`(op(\$(CommonTools.symb_of_symb_list exprs))) >> in
(* If the solution is finite (no infinite sum involved), it goes to *)
(* ser:-finite_part, otherwise to ser:-closed_form. *)
(* In that case, we also update the recurrences and initial values. *)
if <:bool< \$(expr) = 0 >> then
let rec_inits =
<<
proc(xx, u, i)
local poly, deg, ord, recsys, inits, new_inits, eqns, sol;
poly := subs(xx = xx^\$(ram), \$(sol):-poly);
poly := combine(poly, 'power', 'symbolic');
deg := degree(poly, xx);
recsys := \$(ser):-recurrences;
inits := \$(ser):-initial_conditions;
ord := -min(map(a -> op(a) - \$(n),
indets(recsys, specfunc(anything, u))));
eqns := [seq(op(subs(\$(n) = ord + j, recsys)), j = 0 .. deg)];
eqns := subs(
[seq(u[i](j) = coeff(poly, xx, j), j = 0 .. ord + deg)], eqns);
new_inits := remove(a -> evalb(eval(op(2, a), u = 0) = 0), eqns);
eqns := subs(inits, eqns);
sol := solve(eqns, indets(eqns, specfunc(integer, u)));
if [sol] = [] then error("Error 666: no solution found.") end if;
new_inits := select(has, [op(inits), op(sol)],
indets(map2(op, 1, new_inits), specfunc(integer, u)));
recsys := select(has, eval(recsys, u[i] = 0), u);
inits := sort(remove(has, [op(inits), op(new_inits)], u[i]));
return [recsys, inits];
end proc(\$(x_loc), \$(u), \$(i))
>> in
<:unit<
\$(ser):-finite_part := \$(ser):-finite_part + \$(ln) * \$(sol):-poly;
\$(ser):-log_powers := remove(`=`, \$(ser):-log_powers, \$(i));
\$(ser):-recurrences := \$(rec_inits)[1];
\$(ser):-initial_conditions := \$(rec_inits)[2];
>>
else
<:unit<
\$(ser):-closed_form :=
\$(ser):-closed_form + \$(ln) * (\$(sol):-poly + \$(expr));
>>
else
(* If no closed form could be found, a symbolic expression in terms *)
(* of the unknown coefficient sequences u[i](n) is stored in *)
(* ser:-closed_form. *)
let x_part = << \$(x_loc) ^ (\$(n) / \$(ram) + \$(x_expon)) >> in
let smnd = << __DynaMoW_times([\$(u)[\$(i)](\$(n)), \$(x_part)]) >> in
<:unit<
\$(ser):-closed_form :=
\$(ser):-closed_form + \$(ln) * Sum(\$(smnd), \$(n) = 0 .. infinity);
>>
in

(* If this condition is true, it means that the closed form was never *)
(* computed before and that we actually have to do something. *)
if <:bool< \$(ser):-closed_form = undefined >>
then begin
<:unit< \$(ser):-closed_form := 0 >>;
let lps = << \$(ser):-log_powers >> in
let closed_forms =
call_generalTerms
<< \$(ser):-recurrences >> << \$(ser):-initial_conditions >>
u n << \$(x_loc)^(1 / \$(ram)) >> lps in
(* Update the transseries record for all indices i. *)
List.iter
(fun i ->
update_index << \$(lps)[\$(int: i)] >> << \$(closed_forms)[\$(int: i)] >>)
(CommonTools.range <:int< nops(\$(lps)) >>);
end;

(* Now in any case the closed form is available in the record. *)
<:unit<
\$(ser):-finite_part :=
DDMF:-replace_and_combine(\$(ser):-finite_part, \$(x_loc), \$(x));
\$(ser):-closed_form :=
DDMF:-replace_and_combine(\$(ser):-closed_form, \$(x_loc), \$(x));
>>;
let cf = << \$(ser):-finite_part + \$(ser):-closed_form >> in
let recsys, inits =
if <:bool< has(\$(ser):-closed_form, \$(u)) >>
then << \$(ser):-recurrences >>, << \$(ser):-initial_conditions >>
else << [] >>, << [] >> in
let exp_factor = << combine(\$(ser):-exp_factor, 'power', 'symbolic') >> in
<< [\$(exp_factor) * \$(pow_factor) * \$(cf), \$(recsys), \$(inits)] >>

let closed_form_of_linear_comb linear_comb =
let sers = CommonTools.symb_list_of_symb << map(a -> a[2], \$(linear_comb)) >>
and coeffs = << map(a -> a[1], \$(linear_comb)) >> in
let closed_forms =
CommonTools.symb_of_symb_list (List.map closed_form_of_transseries sers) in
let recsys = << map(a -> op(a[2]), \$(closed_forms)) >>
and inits = << map(a -> op(a[3]), \$(closed_forms)) >>
and closed_forms = << map(a -> a[1], \$(closed_forms)) >> in
(* Assemble the linear combination of closed forms. *)
let cf =
<<
`+`(op(zip((a, b) -> `if`(op(0, b) = `+`, map(c -> a * c, b), a * b),
\$(coeffs), \$(closed_forms))))
>> in
let cf = << `if`(op(0, \$(cf)) = `+`, [op(\$(cf))], [\$(cf)]) >> in
(* Sort the terms according to asymptotic behaviour. *)
(* This is not perfectly done since for the sums it's difficult to see. *)
let var = if sers = [] then << x >> else << \$(List.hd sers):-variable >> in
let cf = << selectremove(has, sort(\$(cf), length), Sum) >> in
let cf =
<<
[op(sort(\$(cf)[2],
proc(a, b) local da, db;
da := DDMF:-my_degree(a, ln(\$(var)));
db := DDMF:-my_degree(b, ln(\$(var)));
`if`(da = db, [a, b] = sort([a, b]), da > db)
end proc)), op(\$(cf)[1])]
>> in
(* Simplify coefficients, in particular in front of sums. *)
let cf =
<<
map(a -> `if`(op(0, a) = `*`, ((b, c) -> combine(`*`(op(c))) *
(`*`(op(b))))(selectremove(has, [op(a)], Sum)), a), \$(cf))
>> in
(* Rewrite the parts of cf in terms of __DynaMoW_plus, _times and _frac. *)
let cf =
<<
map(a -> `if`(op(0, a) = `*`, DDMF:-dynamow_times_to_frac(__DynaMoW_times(
map(b -> op(sort(b, (a1, a2) -> `if`(has(a1, Sum) = has(a2, Sum),
[a1, a2] = sort([a1, a2]), has(a2, Sum)))),
[selectremove(b -> not(has(b, indets(\$(var)))), [op(a)])])),
[Sum, \$(var), 1/\$(var)]), a), \$(cf))
>> in
let cf = << `if`(nops(\$(cf)) > 1, __DynaMoW_plus(\$(cf)), `+`(op(\$(cf)))) >> in
(* If there are unknown sequences (no closed form found), rename them, *)
(* as to have consecutive indexing. *)
let seq_names = << map(a -> a[2]:-sequence_name, \$(linear_comb)) >> in
let seq_names =
<< map(a -> `if`(op(0, a) = 'symbol', a, op(0, a)), \$(seq_names)) >> in
if <:bool< nops({op(\$(seq_names))}) = 1 >> then
let rename =
<<
DDMF:-rename_sequences(
[\$(cf), \$(recsys), \$(inits)], \$(seq_names)[1])[1]
>> in
<< \$(rename)[1] >>, << \$(rename)[2] >>, << \$(rename)[3] >>
else cf, recsys, inits

let_service ClosedForm
(linear_comb : any maple)
(notation : string) :

DC.sec_entities * (any maple * any maple) with { title = title } =

let linear_comb = << DDMF:-decode_linear_combination(\$(linear_comb)) >> in
let cf, recsys, inits = closed_form_of_linear_comb linear_comb in
let cf = << DDMF:-maple_to_dynamow(\$(cf)) >> in

<:par<Taylor coefficients: <:dmath< \$(str: notation) = \$(symb: cf). >>>>,
(cf, << DDMF:-encode_linear_combination(\$(linear_comb)) >>)
```

Generated by GNU Enscript 1.6.5.90.