{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Decomposition algorithm" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "3-element Array{DynamicPolynomials.PolyVar{true},1}:\n", " x1\n", " x2\n", " x3" ] }, "execution_count": 1, "metadata": {}, "output_type": "execute_result" } ], "source": [ "using TensorDec\n", "using LinearAlgebra\n", "X = @ring x1 x2 x3" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We want to find a sparse representation of the following series known up to degree 3:" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "4.0dx1*dx3 + 15.0dx2 + 6.0dx3 + 20.0dx1*dx2 + 6.0dx3^3 + 43.0dx2^2dx3 - 26.0dx1^3 + 129.0dx2^3 + 30.0dx1^2dx2 + 15.0dx2*dx3 + 20.0dx1*dx2*dx3 + 6.0dx1^2 + 6.0dx3^2 + 4.0dx1 + 43.0dx2^2 + 6.0dx1^2dx3 + 4.0dx1*dx3^2 + 72.0dx1*dx2^2 + 15.0dx2*dx3^2 + 6.0" ] }, "execution_count": 2, "metadata": {}, "output_type": "execute_result" } ], "source": [ "sigma = dual(6.0 + 4.0*x1 + 15.0*x2 + 6.0*x3 + 6.0*x1^2 + 20.0*x1*x2 + 4.0*x1*x3 + 43.0*x2^2 + 15.0*x2*x3 + 6.0*x3^2 - 26.0*x1^3 + 30.0*x1^2*x2 + 6.0*x1^2*x3 + 72.0*x1*x2^2 + 20.0*x1*x2*x3 + 4.0*x1*x3^2 + 129.0*x2^3 + 43.0*x2^2*x3 + 15.0*x2*x3^2 + 6.0*x3^3)" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "10-element Array{DynamicPolynomials.Monomial{true},1}:\n", " 1 \n", " x1 \n", " x2 \n", " x3 \n", " x1² \n", " x1x2\n", " x1x3\n", " x2² \n", " x2x3\n", " x3² " ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "L1 = monoms(X,1)\n", "L2 = monoms(X,2)" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "4×10 Array{Float64,2}:\n", " 6.0 4.0 15.0 6.0 6.0 20.0 4.0 43.0 15.0 6.0\n", " 4.0 6.0 20.0 4.0 -26.0 30.0 6.0 72.0 20.0 4.0\n", " 15.0 20.0 43.0 15.0 30.0 72.0 20.0 129.0 43.0 15.0\n", " 6.0 4.0 15.0 6.0 6.0 20.0 4.0 43.0 15.0 6.0" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "H = hankel(sigma,L1,L2)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The rank of $H_{\\sigma}$ will give us an idea on the dimension of $\\mathcal{A}_\\sigma$." ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "3" ] }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" } ], "source": [ "rank(H)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We check that $\\{1, x_1, x_2\\}$ is a basis of $\\mathcal{A}_\\sigma$: " ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "3-element Array{DynamicPolynomials.Monomial{true},1}:\n", " 1 \n", " x1\n", " x2" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "B0 = L1[1:3]" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "3×3 Array{Float64,2}:\n", " 6.0 4.0 15.0\n", " 4.0 6.0 20.0\n", " 15.0 20.0 43.0" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "H0 = hankel(sigma, B0, B0)" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "3" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "rank(H0)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let us compute the shifted (truncated) Hankel operators." ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "3×3 Array{Float64,2}:\n", " 4.0 6.0 20.0\n", " 6.0 -26.0 30.0\n", " 20.0 30.0 72.0" ] }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" } ], "source": [ "H1 = hankel(sigma, B0, B0*x1)\n", "H2 = hankel(sigma, B0, B0*x2)\n", "H3 = hankel(sigma, B0, B0*x3);\n", "H = [H1,H2,H3]\n", "H[1]" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "3×3 Array{Float64,2}:\n", " 1.11022e-16 9.14286 -0.571429\n", " 1.0 3.85714 1.57143 \n", " -1.11022e-16 -4.28571 1.14286 " ] }, "execution_count": 10, "metadata": {}, "output_type": "execute_result" } ], "source": [ "M = [ H0^(-1)*H[i] for i in 1:3 ]\n", "M[1]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The eigenvalues and eigenvectors of $M_{x_1}$ are" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We deduce the operators of multiplication by the variables in the basis $B_0$:" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "Eigen{Float64,Float64,Array{Float64,2},Array{Float64,1}}\n", "eigenvalues:\n", "3-element Array{Float64,1}:\n", " -0.9999999999999991\n", " 4.000000000000002 \n", " 2.000000000000002 \n", "eigenvectors:\n", "3×3 Array{Float64,2}:\n", " 0.963087 -0.811107 -0.762001\n", " -0.120386 -0.324443 -0.127 \n", " -0.240772 0.486664 0.635001" ] }, "execution_count": 11, "metadata": {}, "output_type": "execute_result" } ], "source": [ "v, E = eigen(M[1])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The matrices $M_{x_i}$ are diagonal in this basis:" ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "3×3 Array{Float64,2}:\n", " -1.0 -6.99441e-15 -3.66374e-15\n", " 4.21885e-15 4.0 -4.44089e-15\n", " -4.66294e-15 -3.9968e-15 2.0 " ] }, "execution_count": 12, "metadata": {}, "output_type": "execute_result" } ], "source": [ "D = [E^(-1)*M[i]*E for i in 1:3]\n", "D[1]" ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "3×3 Array{Float64,2}:\n", " 1.0 -4.44089e-16 -1.44329e-15\n", " 8.88178e-16 2.0 2.66454e-15\n", " -3.55271e-15 2.66454e-15 3.0 " ] }, "execution_count": 13, "metadata": {}, "output_type": "execute_result" } ], "source": [ "D[2]" ] }, { "cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "3×3 Array{Float64,2}:\n", " 1.0 3.33067e-16 1.11022e-16\n", " -9.4369e-16 1.0 6.66134e-16\n", " 5.55112e-16 -6.66134e-16 1.0 " ] }, "execution_count": 14, "metadata": {}, "output_type": "execute_result" } ], "source": [ "D[3]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Looking at the corresponding terms on the diagonal, we get the coordinates of the points $\\Xi$:" ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "3×3 Array{Float64,2}:\n", " -1.0 4.0 2.0\n", " 1.0 2.0 3.0\n", " 1.0 1.0 1.0" ] }, "execution_count": 15, "metadata": {}, "output_type": "execute_result" } ], "source": [ "Xi = [ D[i][j,j] for i in 1:3, j in 1:3]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We normalize the eigenvectors by $v_i \\over v_i(\\xi_i)$ and get the interpolation polynomials at the points $\\xi_i$:" ] }, { "cell_type": "code", "execution_count": 16, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "3-element Array{DynamicPolynomials.Polynomial{true,Float64},1}:\n", " -0.14285714285714324x1 - 0.2857142857142862x2 + 1.142857142857143 \n", " 0.28571428571428614x1 - 0.4285714285714279x2 + 0.7142857142857121 \n", " -0.14285714285714332x1 + 0.7142857142857134x2 - 0.8571428571428543" ] }, "execution_count": 16, "metadata": {}, "output_type": "execute_result" } ], "source": [ "Dg = E'*vcat(fill(1.,1,3), Xi[1:2,:])\n", "E = E*Dg^(-1)\n", "U = E'*B0" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We deduce the weights $w_i=\\sigma(u_i)$:" ] }, { "cell_type": "code", "execution_count": 17, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "3×1 Array{Float64,2}:\n", " 1.999999999999992 \n", " -1.0000000000000018\n", " 5.000000000000002 " ] }, "execution_count": 17, "metadata": {}, "output_type": "execute_result" } ], "source": [ "w = hankel(sigma, U, [L1[1]])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Using the command decompose, we can get directly the same decomposition: " ] }, { "cell_type": "code", "execution_count": 18, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "([-1.0, 5.0, 2.0], [4.0 2.0 -1.0; 2.0 3.0 1.0; 1.0 1.0 1.0])" ] }, "execution_count": 18, "metadata": {}, "output_type": "execute_result" } ], "source": [ "w, Xi = decompose(sigma)" ] }, { "cell_type": "code", "execution_count": 19, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "3×3 Array{Float64,2}:\n", " 4.0 2.0 -1.0\n", " 2.0 3.0 1.0\n", " 1.0 1.0 1.0" ] }, "execution_count": 19, "metadata": {}, "output_type": "execute_result" } ], "source": [ "Xi" ] }, { "cell_type": "code", "execution_count": 20, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "3-element Array{Float64,1}:\n", " -1.0000000000000129\n", " 5.000000000000011 \n", " 1.9999999999999998" ] }, "execution_count": 20, "metadata": {}, "output_type": "execute_result" } ], "source": [ "w" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The series decomposes as $2 \\mathfrak{e}_{(-1,1,1)} + 5 \\mathfrak{e}_{(2,3,1)} - \\mathfrak{e}_{(4,2,1)}$." ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "anaconda-cloud": {}, "kernelspec": { "display_name": "Julia 1.0.0", "language": "julia", "name": "julia-1.0" }, "language_info": { "file_extension": ".jl", "mimetype": "application/julia", "name": "julia", "version": "1.0.0" } }, "nbformat": 4, "nbformat_minor": 1 }