Documentation‎ > ‎NCSOStoolsdemo‎ > ‎

### NCSOStoolsdemo - Rationalization

 %************************************************************************** % Rationalization                                                         * %************************************************************************** % In this section we will look at the module RprojRldlt which tries to find % an exact rational positive semidefinite solution of the problem A*X=b % from the floating point solution XX. The algorithm tries to find an % appropiate rational approximation of the matrix XX which is after the % projection on the hyperplane A*X=b a positive semidefinite matrix. % Module has also an option to try to compute an exact rational LDU % decomposition of a rational positive semidefinite solution of the problem % A*XX=b. Otherwise positive semidefiniteness relies on numerical % computation of eigenvalues. % This is usefull whenever we want to find a rational Gram matrix an % therefore a rational sum of hermitian squares decomposition, i.e. in % module NCsos, NCmin, NCcycSos, NCcycMin ... % ------------------------------------------------------------------------- % First we construct some symbolic noncommuting variables and polynomial % ------------------------------------------------------------------------- >> NCvars x y >> f = 1 + 3*x^2*y*x^2*y + 7*x^2*y^2*x^2 + 2*x*y*x + 6*x*y*x*y*x^2; % ------------------------------------------------------------------------- % With NCcycSos we can get numerical solution and therefore numerical % (inexact) sum of hermitian squares decomposition. This will be our % starting point in finding exact solution. % For the rationalization is in most cases better to take objective % function for the SDP solver to be zero; so we set this in parameter. % ------------------------------------------------------------------------- >> params.obj = 0; [IsCycEq,XX,base,sohs,g,SDP_data] = NCcycSos(f, params) ***** NCSOStools: module NCcycSos started ***** Input polynomial has (max) degree 6 and min degree 0. Detected 5 monomials in 2 variables. There are 127 monomials in 2 variables of degree at most 6. Preprocessing the input polynomial ... Using alpha degree to construct the monomial vector. Keeping 4 monomials in the monomial vector. Computing cyclically equivalent products ... done. Preparing linear constraints ... Number of linear constraints: 5. Starting SDP solver ... ...Residual norm: 4.6902e-014 ***** ALERT: SDP solver ran into minor numerical problems.                        ***** ***** Nevertheless, according to feasratio (0.8878) given solution MIGHT be ok.   ***** Computing SOHS decomposition ... done. Found SOHS decomposition with 4 factors. WARNING! SDP solver returned some numerical problems. Check messages! IsCycEq = 1 XX = 1.0000    0.0439    0.5087    0.4474      0.0439    2.1604    0.6947   -1.7065      0.5087    0.6947    6.4130    2.3053      0.4474   -1.7065    2.3053    4.8396 base = ''        'x*x*y'        'x*y*x'        'y*x*x' sohs = 1+0.0439*x^2*y+0.5087*x*y*x+0.4474*y*x^2            1.469*x^2*y+0.4577*x*y*x-1.175*y*x^2                         2.438*x*y*x+1.073*y*x^2                                     1.452*y*x^2 g = 1+0.4913*x^2*y - 1.7065222*x^2*y*x^2*y + 2.3052978*x^2*y*x*y*x + 4.8395461*x^2*y^2*x^2 + 1.0174*x*y*x + 0.69478477*x*y*x^3*y + 6.4130842*x*y*x^2*y*x + 2.3052978*x*y*x*y*x^2 + 0.4913*y*x^2 + 2.1604759*y*x^4*y + 0.69478477*y*x^3*y*x - 1.7065222*y*x^2*y*x^2 SDP_data =        A: [5x16 double]        b: [5x1 double]        C: [16x1 double]        K: [1x1 struct]     pars: [1x1 struct] % ------------------------------------------------------------------------- % Now we can use the algorithm in RprojRldlt to try to find an appropiate % rational approximation of the matrix XX which is after the projection on % the hyperplane A*X=b a positive semidefinite matrix; if we succeed in % finding rational Gram matrix we can therefore get a rational sum of % hermitian squares decomposition. The next step depends heavily on the % 'quality' of the Gram matrix XX returned by the SDP solver. Here we % measure 'quality' as large as possible minimal eigenvalue of XX. % ------------------------------------------------------------------------- >> Xrat = RprojRldlt(XX,SDP_data,true); Xrat{1}, Xrat{2} ***** NCSOStools: module RProjRldlt started ***** delta (min eig X): 7.373264e-001 epsilon (residual norm): 4.690190e-014  i :    1/2^i      tau_i      min_eig  0 :  1.00E+000  8.79E-001  4.49E-001 Rational projection and LDLT decomposition found. ***** NCSOStools: module RProjRldlt completed successfully ***** ans = 1     0     1     0       0     2     1    -5       1     1    19     2       0    -5     2     5 ans = 1     1     1     1       1     1     1     3       1     1     3     1       1     3     1     1 % ------------------------------------------------------------------------- % Now we got an exact rational Gram matrix Xrat where Xrat{1} are % numerators and Xrat{2} adherent denominators. If we also want a rational % sum of hermitian squares decomposition we need to slightly change a call % to the module RprojRldlt to get a rational LDU decomposition of Xrat, % more precise we get a rational decomposition of the form % Xrat=(P*L)*D*(P*L)^T. % ------------------------------------------------------------------------- >> [Xrat,L,D,P,err] = RprojRldlt(XX,SDP_data,true); ***** NCSOStools: module RProjRldlt started ***** delta (min eig X): 7.373264e-001 epsilon (residual norm): 4.690190e-014  i :    1/2^i      tau_i      min_eig  0 :  1.00E+000  8.79E-001  4.49E-001 Rational projection and LDLT decomposition found. ***** NCSOStools: module RProjRldlt completed successfully ***** >> base base = ''        'x*x*y'        'x*y*x'        'y*x*x' >> P P = 0     0     0     1     0     0     1     0     1     0     0     0     0     1     0     0 >> D{1}, D{2} ans = 19     0     0     0        0    83     0     0        0     0   704     0        0     0     0   509 ans = 3     1     1     1       1    19     1     1       1     1   747     1       1     1     1   704 >> L{1}, L{2} ans = 1     0     0     0       6     1     0     0       3  -113     1     0       3    -6  -225     1 ans =  1     1     1     1       19     1     1     1       19   249     1     1       19    83   704     1 % If RprojRldlt was succesfull we gain an exact rational in the following % way: % - first we permute a base vector according to the permutation matrix P % - each column of the matrix L corresponds to the coefficients of the one % hermitian square in the decomposition multiplied by the corresponding % element in the matrix D.