Computation of Scalar Products

By now, we are able to construct orthogonal polynomials, and to construct quadrature rules for a given nonnegative weight function, respectively. Now we combine both ideas to solve integrals involving the orthogonal polynomials

\[\langle \phi_{i_1} \phi_{i_2} \cdots \phi_{i_{m-1}}, \phi_{i_m} \rangle = \int \phi_{i_1}(t) \phi_{i_2}(t) \cdots \phi_{i_{m-1}}(t) \phi_{i_m}(t) w(t) \mathrm{d} t,\]

both for the univariate and multivariate case. The integrand is a polynomial (possibly multivariate) that can be solved exactly with the appropriate Gauss quadrature rules.

Note

To simplify notation we drop the integration interval. It is clear from the context.

Univariate Polynomials

Classical Polynomials

Let's begin with a univariate basis for some classical orthogonal polynomial

using PolyChaos
deg, n = 4, 20
s_α, s_β = 2.1, 3.2
opq = Beta01OrthoPoly(deg, s_α, s_β; Nrec=n, addQuadrature=true)
Beta01OrthoPoly{Vector{Float64}, Beta01Measure, Quad{Float64, Vector{Float64}}}(4, [0.39622641509433965, 0.45308865339881105, 0.4732655766681396, 0.482729089351984, 0.4879233481934926, 0.49108064278342917, 0.49314292190864784, 0.4945640779897571, 0.4955849084142159, 0.4963428640512603, 0.49692106736331404, 0.4973721930243669, 0.49773093800397555, 0.4980209139779084, 0.4982586420286101, 0.4984559630451982, 0.4986215434156863, 0.4987618443402984, 0.4988817625639983, 0.4989850639437675], [1.0, 0.03797318144060756, 0.049528193464418745, 0.054496397032519545, 0.05707494398498378, 0.05858224318366976, 0.05953877522858778, 0.06018346810103027, 0.0606384719605267, 0.06097150420824787, 0.061222558056958205, 0.06141648535528466, 0.061569388295544744, 0.061692071403646374, 0.06179200321800959, 0.061874480034347076, 0.06194334217794196, 0.06200142915302492, 0.06205087705273863, 0.062093317746215744], Beta01Measure(PolyChaos.var"#102#103"{Float64, Float64}(2.1, 3.2), (0.0, 1.0), false, 2.1, 3.2), Quad{Float64, Vector{Float64}}("golubwelsch", 19, [0.00878890525374891, 0.028450106116514506, 0.05852125536671723, 0.09833813997192814, 0.14702350691009355, 0.20350504892083368, 0.266538831062297, 0.3347366435399401, 0.40659656306145764, 0.48053602819819125, 0.5549266957500689, 0.6281303108338084, 0.6985348047519305, 0.7645898365922597, 0.8248410256009169, 0.8779622119533655, 0.9227853738738273, 0.9583291958712525, 0.9838388819882817], [0.00115167784765143, 0.007007259624729887, 0.020314892689531894, 0.04140255942833847, 0.06777171117046146, 0.09471990538153915, 0.11673118518989155, 0.12913645448496422, 0.12947566141190417, 0.11812235962361355, 0.09800772591714006, 0.07359961541684255, 0.04953526816082873, 0.029389573634464723, 0.014963964676256315, 0.006252710931950127, 0.0019796669915308576, 0.00040351473743713024, 3.429268092282879e-5]))

By setting addQuadrature = true (which is default), an $n$-point Gauss quadrature rule is create relative to the underlying measure opq.measure, where $n$ is the number of recurrence coefficients stored in opq.α and opq.β.

To compute the squared norms

\[\| \phi_k \|^2 = \langle \phi_k, \phi_k \rangle = \int \phi_k(t) \phi_k(t) w(t) \mathrm{d} t\]

of the basis we call computeSP2()

normsq = computeSP2(opq)
5-element Vector{Float64}:
 1.0
 0.03797318144060756
 0.0018807430768498865
 0.00010249372143217383
 5.849823409553853e-6

For the general case

\[\langle \phi_{i_1} \phi_{i_2} \cdots \phi_{i_{m-1}}, \phi_{i_m} \rangle = \int \phi_{i_1}(t) \phi_{i_2}(t) \cdots \phi_{i_{m-1}}(t) \phi_{i_m}(t) w(t) \mathrm{d} t,\]

there exists a type Tensor that requires only two arguments: the dimension $m \geq 1$, and an AbstractOrthoPoly

m = 3
t = Tensor(3,opq)
Tensor{Beta01OrthoPoly{Vector{Float64}, Beta01Measure, Quad{Float64, Vector{Float64}}}}(3,   [1 ]  =  1.0
  [6 ]  =  0.0379732
  [11]  =  0.00188074
  [16]  =  0.000102494
  [21]  =  5.84982e-6
  [22]  =  0.00215924
  [23]  =  0.00188074
        ⋮
  [45]  =  5.84982e-6
  [48]  =  7.80614e-6
  [49]  =  7.09802e-7
  [53]  =  4.73123e-7
  [64]  =  9.40423e-7
  [65]  =  5.08904e-7
  [69]  =  6.53232e-8
  [85]  =  3.55404e-8, PolyChaos.var"#getfun#45"{Int64, Beta01OrthoPoly{Vector{Float64}, Beta01Measure, Quad{Float64, Vector{Float64}}}, SparseArrays.SparseVector{Float64, Int64}}(3, Beta01OrthoPoly{Vector{Float64}, Beta01Measure, Quad{Float64, Vector{Float64}}}(4, [0.39622641509433965, 0.45308865339881105, 0.4732655766681396, 0.482729089351984, 0.4879233481934926, 0.49108064278342917, 0.49314292190864784, 0.4945640779897571, 0.4955849084142159, 0.4963428640512603, 0.49692106736331404, 0.4973721930243669, 0.49773093800397555, 0.4980209139779084, 0.4982586420286101, 0.4984559630451982, 0.4986215434156863, 0.4987618443402984, 0.4988817625639983, 0.4989850639437675], [1.0, 0.03797318144060756, 0.049528193464418745, 0.054496397032519545, 0.05707494398498378, 0.05858224318366976, 0.05953877522858778, 0.06018346810103027, 0.0606384719605267, 0.06097150420824787, 0.061222558056958205, 0.06141648535528466, 0.061569388295544744, 0.061692071403646374, 0.06179200321800959, 0.061874480034347076, 0.06194334217794196, 0.06200142915302492, 0.06205087705273863, 0.062093317746215744], Beta01Measure(PolyChaos.var"#102#103"{Float64, Float64}(2.1, 3.2), (0.0, 1.0), false, 2.1, 3.2), Quad{Float64, Vector{Float64}}("golubwelsch", 19, [0.00878890525374891, 0.028450106116514506, 0.05852125536671723, 0.09833813997192814, 0.14702350691009355, 0.20350504892083368, 0.266538831062297, 0.3347366435399401, 0.40659656306145764, 0.48053602819819125, 0.5549266957500689, 0.6281303108338084, 0.6985348047519305, 0.7645898365922597, 0.8248410256009169, 0.8779622119533655, 0.9227853738738273, 0.9583291958712525, 0.9838388819882817], [0.00115167784765143, 0.007007259624729887, 0.020314892689531894, 0.04140255942833847, 0.06777171117046146, 0.09471990538153915, 0.11673118518989155, 0.12913645448496422, 0.12947566141190417, 0.11812235962361355, 0.09800772591714006, 0.07359961541684255, 0.04953526816082873, 0.029389573634464723, 0.014963964676256315, 0.006252710931950127, 0.0019796669915308576, 0.00040351473743713024, 3.429268092282879e-5])),   [1 ]  =  1.0
  [6 ]  =  0.0379732
  [11]  =  0.00188074
  [16]  =  0.000102494
  [21]  =  5.84982e-6
  [22]  =  0.00215924
  [23]  =  0.00188074
        ⋮
  [45]  =  5.84982e-6
  [48]  =  7.80614e-6
  [49]  =  7.09802e-7
  [53]  =  4.73123e-7
  [64]  =  9.40423e-7
  [65]  =  5.08904e-7
  [69]  =  6.53232e-8
  [85]  =  3.55404e-8), Beta01OrthoPoly{Vector{Float64}, Beta01Measure, Quad{Float64, Vector{Float64}}}(4, [0.39622641509433965, 0.45308865339881105, 0.4732655766681396, 0.482729089351984, 0.4879233481934926, 0.49108064278342917, 0.49314292190864784, 0.4945640779897571, 0.4955849084142159, 0.4963428640512603, 0.49692106736331404, 0.4973721930243669, 0.49773093800397555, 0.4980209139779084, 0.4982586420286101, 0.4984559630451982, 0.4986215434156863, 0.4987618443402984, 0.4988817625639983, 0.4989850639437675], [1.0, 0.03797318144060756, 0.049528193464418745, 0.054496397032519545, 0.05707494398498378, 0.05858224318366976, 0.05953877522858778, 0.06018346810103027, 0.0606384719605267, 0.06097150420824787, 0.061222558056958205, 0.06141648535528466, 0.061569388295544744, 0.061692071403646374, 0.06179200321800959, 0.061874480034347076, 0.06194334217794196, 0.06200142915302492, 0.06205087705273863, 0.062093317746215744], Beta01Measure(PolyChaos.var"#102#103"{Float64, Float64}(2.1, 3.2), (0.0, 1.0), false, 2.1, 3.2), Quad{Float64, Vector{Float64}}("golubwelsch", 19, [0.00878890525374891, 0.028450106116514506, 0.05852125536671723, 0.09833813997192814, 0.14702350691009355, 0.20350504892083368, 0.266538831062297, 0.3347366435399401, 0.40659656306145764, 0.48053602819819125, 0.5549266957500689, 0.6281303108338084, 0.6985348047519305, 0.7645898365922597, 0.8248410256009169, 0.8779622119533655, 0.9227853738738273, 0.9583291958712525, 0.9838388819882817], [0.00115167784765143, 0.007007259624729887, 0.020314892689531894, 0.04140255942833847, 0.06777171117046146, 0.09471990538153915, 0.11673118518989155, 0.12913645448496422, 0.12947566141190417, 0.11812235962361355, 0.09800772591714006, 0.07359961541684255, 0.04953526816082873, 0.029389573634464723, 0.014963964676256315, 0.006252710931950127, 0.0019796669915308576, 0.00040351473743713024, 3.429268092282879e-5])))

To get the desired entries, Tensor comes with a get() function that is called for some index $a \in \mathbb{N}_0^m$ that has the entries $a = [i_1, i_2, \dots, i_m]$. For example

t.get([1,2,3])
0.00010249372143217402

Or using comprehension

T = [ t.get([i1,i2,i3]) for i1=0:dim(opq)-1, i2=0:dim(opq)-1, i3=0:dim(opq)-1]
5×5×5 Array{Float64, 3}:
[:, :, 1] =
 1.0  0.0        0.0         0.0          0.0
 0.0  0.0379732  0.0         0.0          0.0
 0.0  0.0        0.00188074  0.0          0.0
 0.0  0.0        0.0         0.000102494  0.0
 0.0  0.0        0.0         0.0          5.84982e-6

[:, :, 2] =
 0.0        0.0379732   0.0          0.0          0.0
 0.0379732  0.00215924  0.00188074   0.0          0.0
 0.0        0.00188074  0.000144891  0.000102494  0.0
 0.0        0.0         0.000102494  8.86598e-6   5.84982e-6
 0.0        0.0         0.0          5.84982e-6   5.36411e-7

[:, :, 3] =
 0.0         0.0          0.00188074   0.0          0.0
 0.0         0.00188074   0.000144891  0.000102494  0.0
 0.00188074  0.000144891  0.000127149  1.0934e-5    5.84982e-6
 0.0         0.000102494  1.0934e-5    7.80614e-6   7.09802e-7
 0.0         0.0          5.84982e-6   7.09802e-7   4.73123e-7

[:, :, 4] =
 0.0          0.0          0.0          0.000102494  0.0
 0.0          0.0          0.000102494  8.86598e-6   5.84982e-6
 0.0          0.000102494  1.0934e-5    7.80614e-6   7.09802e-7
 0.000102494  8.86598e-6   7.80614e-6   9.40423e-7   5.08904e-7
 0.0          5.84982e-6   7.09802e-7   5.08904e-7   6.53232e-8

[:, :, 5] =
 0.0         0.0         0.0         0.0         5.84982e-6
 0.0         0.0         0.0         5.84982e-6  5.36411e-7
 0.0         0.0         5.84982e-6  7.09802e-7  4.73123e-7
 0.0         5.84982e-6  7.09802e-7  5.08904e-7  6.53232e-8
 5.84982e-6  5.36411e-7  4.73123e-7  6.53232e-8  3.55404e-8

Notice that we can cross-check the results.

using LinearAlgebra
normsq == diag(T[:,:,1]) == diag(T[:,1,:]) == diag(T[1,:,:])
true

Also, normsq can be computed analogously in Tensor format

t2 = Tensor(2, opq)
normsq == [ t2.get([i, i]) for i in 0:dim(opq)-1]
true

Arbitrary Weights

Of course, the type OrthoPoly can be constructed for arbitrary weights $w(t)$. In this case we have to compute the orthogonal basis and the respective quadrature rule. Let's re-work the above example by hand.

using SpecialFunctions
supp = (0, 1)
w(t) = (t^(s_α-1)*(1-t)^(s_β-1)/SpecialFunctions.beta(s_α,s_β))
my_meas = Measure("my_meas", w, supp, false)
my_opq = OrthoPoly("my_op", deg, my_meas; Nrec=n, addQuadrature = true)
OrthoPoly{Vector{Float64}, Measure, Quad{Float64, Vector{Float64}}}("my_op", 4, [0.39622641509448137, 0.45308865339922233, 0.47326557666900376, 0.4827290893535496, 0.48792334819610145, 0.4910806427875327, 0.4931429219148361, 0.4945640779987818, 0.495584908427019, 0.4963428640690031, 0.49692106738740605, 0.4973721930564992, 0.4977309380461525, 0.4980209140324845, 0.49825864209832377, 0.49845596313321483, 0.4986215435256319, 0.4987618444763114, 0.4988817627307632, 0.4989850641465735], [0.9999999999996537, 0.037973181440564815, 0.049528193464296905, 0.054496397032274234, 0.05707494398455889, 0.05858224318299423, 0.059538775227572845, 0.06018346809956606, 0.060638471958479156, 0.060971504205455586, 0.0612225580532287, 0.061416485350391216, 0.06156938828922262, 0.06169207139558928, 0.06179200320786561, 0.06187448002171521, 0.06194334216236786, 0.062001429133996595, 0.06205087702968264, 0.06209331771849251], Measure("my_meas", Main.w, (0.0, 1.0), false, Dict{Any, Any}()), Quad{Float64, Vector{Float64}}("golubwelsch", 19, [0.008788905307716177, 0.02845010617907043, 0.0585212554299101, 0.09833814003385198, 0.14702350696982944, 0.20350504897777927, 0.2665388311160098, 0.33473664359009403, 0.40659656310782666, 0.4805360282406439, 0.5549266957885687, 0.6281303108684106, 0.6985348047827785, 0.7645898366195809, 0.824841025625018, 0.877962211974623, 0.9227853738926778, 0.9583291958881769, 0.9838388820037522], [0.001151677857036659, 0.007007259641384841, 0.0203148927103142, 0.0414025594490785, 0.06777171118707008, 0.09471990539100715, 0.11673118519090996, 0.1291364544780523, 0.12947566139904687, 0.11812235960758868, 0.09800772590075595, 0.07359961540231305, 0.04953526814944095, 0.02938957362654602, 0.0149639646714071, 0.00625271092938539, 0.001979666990405483, 0.00040351473706229986, 3.429268084906047e-5]))

Now we can compute the squared norms $\| \phi_k \|^2$

my_normsq = computeSP2(my_opq)
5-element Vector{Float64}:
 0.9999999999996537
 0.037973181440551666
 0.0018807430768424916
 0.00010249372143130948
 5.849823409460971e-6

And the tensor

my_t = Tensor(m, my_opq)
my_T = [ my_t.get([i1,i2,i3]) for i1=0:dim(opq)-1,i2=0:dim(opq)-1,i3=0:dim(opq)-1]
5×5×5 Array{Float64, 3}:
[:, :, 1] =
 1.0  0.0        0.0         0.0          0.0
 0.0  0.0379732  0.0         0.0          0.0
 0.0  0.0        0.00188074  0.0          0.0
 0.0  0.0        0.0         0.000102494  0.0
 0.0  0.0        0.0         0.0          5.84982e-6

[:, :, 2] =
 0.0        0.0379732   0.0          0.0          0.0
 0.0379732  0.00215924  0.00188074   0.0          0.0
 0.0        0.00188074  0.000144891  0.000102494  0.0
 0.0        0.0         0.000102494  8.86598e-6   5.84982e-6
 0.0        0.0         0.0          5.84982e-6   5.36411e-7

[:, :, 3] =
 0.0         0.0          0.00188074   0.0          0.0
 0.0         0.00188074   0.000144891  0.000102494  0.0
 0.00188074  0.000144891  0.000127149  1.0934e-5    5.84982e-6
 0.0         0.000102494  1.0934e-5    7.80614e-6   7.09802e-7
 0.0         0.0          5.84982e-6   7.09802e-7   4.73123e-7

[:, :, 4] =
 0.0          0.0          0.0          0.000102494  0.0
 0.0          0.0          0.000102494  8.86598e-6   5.84982e-6
 0.0          0.000102494  1.0934e-5    7.80614e-6   7.09802e-7
 0.000102494  8.86598e-6   7.80614e-6   9.40423e-7   5.08904e-7
 0.0          5.84982e-6   7.09802e-7   5.08904e-7   6.53232e-8

[:, :, 5] =
 0.0         0.0         0.0         0.0         5.84982e-6
 0.0         0.0         0.0         5.84982e-6  5.36411e-7
 0.0         0.0         5.84982e-6  7.09802e-7  4.73123e-7
 0.0         5.84982e-6  7.09802e-7  5.08904e-7  6.53232e-8
 5.84982e-6  5.36411e-7  4.73123e-7  6.53232e-8  3.55404e-8

Let's compare the results:

abs.(normsq-my_normsq)
5-element Vector{Float64}:
 3.462785613805863e-13
 5.5892790395972725e-14
 7.394909337654632e-15
 8.643530769597563e-16
 9.28813978311703e-17
norm(T-my_T)
3.600949810975133e-13
Note

The possibility to create quadrature rules for arbitrary weights should be reserved to cases different from classical ones.

Multivariate Polynomials

For multivariate polynomials the syntax for Tensor is very much alike, except that we are dealing with the type MultiOrthoPoly now.

mop = MultiOrthoPoly([opq, my_opq], deg)
MultiOrthoPoly{ProductMeasure, Quad{Float64, Vector{Float64}}, Vector{AbstractOrthoPoly{M, Quad{Float64, Vector{Float64}}} where M<:AbstractMeasure}}(["Beta01OrthoPoly{Vector{Float64}, Beta01Measure, Quad{Float64, Vector{Float64}}}", "my_op"], 4, 15, [0 0; 1 0; … ; 1 3; 0 4], ProductMeasure(PolyChaos.var"#w#39"{Vector{AbstractOrthoPoly{M, Quad{Float64, Vector{Float64}}} where M<:AbstractMeasure}}(AbstractOrthoPoly{M, Quad{Float64, Vector{Float64}}} where M<:AbstractMeasure[Beta01OrthoPoly{Vector{Float64}, Beta01Measure, Quad{Float64, Vector{Float64}}}(4, [0.39622641509433965, 0.45308865339881105, 0.4732655766681396, 0.482729089351984, 0.4879233481934926, 0.49108064278342917, 0.49314292190864784, 0.4945640779897571, 0.4955849084142159, 0.4963428640512603, 0.49692106736331404, 0.4973721930243669, 0.49773093800397555, 0.4980209139779084, 0.4982586420286101, 0.4984559630451982, 0.4986215434156863, 0.4987618443402984, 0.4988817625639983, 0.4989850639437675], [1.0, 0.03797318144060756, 0.049528193464418745, 0.054496397032519545, 0.05707494398498378, 0.05858224318366976, 0.05953877522858778, 0.06018346810103027, 0.0606384719605267, 0.06097150420824787, 0.061222558056958205, 0.06141648535528466, 0.061569388295544744, 0.061692071403646374, 0.06179200321800959, 0.061874480034347076, 0.06194334217794196, 0.06200142915302492, 0.06205087705273863, 0.062093317746215744], Beta01Measure(PolyChaos.var"#102#103"{Float64, Float64}(2.1, 3.2), (0.0, 1.0), false, 2.1, 3.2), Quad{Float64, Vector{Float64}}("golubwelsch", 19, [0.00878890525374891, 0.028450106116514506, 0.05852125536671723, 0.09833813997192814, 0.14702350691009355, 0.20350504892083368, 0.266538831062297, 0.3347366435399401, 0.40659656306145764, 0.48053602819819125, 0.5549266957500689, 0.6281303108338084, 0.6985348047519305, 0.7645898365922597, 0.8248410256009169, 0.8779622119533655, 0.9227853738738273, 0.9583291958712525, 0.9838388819882817], [0.00115167784765143, 0.007007259624729887, 0.020314892689531894, 0.04140255942833847, 0.06777171117046146, 0.09471990538153915, 0.11673118518989155, 0.12913645448496422, 0.12947566141190417, 0.11812235962361355, 0.09800772591714006, 0.07359961541684255, 0.04953526816082873, 0.029389573634464723, 0.014963964676256315, 0.006252710931950127, 0.0019796669915308576, 0.00040351473743713024, 3.429268092282879e-5])), OrthoPoly{Vector{Float64}, Measure, Quad{Float64, Vector{Float64}}}("my_op", 4, [0.39622641509448137, 0.45308865339922233, 0.47326557666900376, 0.4827290893535496, 0.48792334819610145, 0.4910806427875327, 0.4931429219148361, 0.4945640779987818, 0.495584908427019, 0.4963428640690031, 0.49692106738740605, 0.4973721930564992, 0.4977309380461525, 0.4980209140324845, 0.49825864209832377, 0.49845596313321483, 0.4986215435256319, 0.4987618444763114, 0.4988817627307632, 0.4989850641465735], [0.9999999999996537, 0.037973181440564815, 0.049528193464296905, 0.054496397032274234, 0.05707494398455889, 0.05858224318299423, 0.059538775227572845, 0.06018346809956606, 0.060638471958479156, 0.060971504205455586, 0.0612225580532287, 0.061416485350391216, 0.06156938828922262, 0.06169207139558928, 0.06179200320786561, 0.06187448002171521, 0.06194334216236786, 0.062001429133996595, 0.06205087702968264, 0.06209331771849251], Measure("my_meas", Main.w, (0.0, 1.0), false, Dict{Any, Any}()), Quad{Float64, Vector{Float64}}("golubwelsch", 19, [0.008788905307716177, 0.02845010617907043, 0.0585212554299101, 0.09833814003385198, 0.14702350696982944, 0.20350504897777927, 0.2665388311160098, 0.33473664359009403, 0.40659656310782666, 0.4805360282406439, 0.5549266957885687, 0.6281303108684106, 0.6985348047827785, 0.7645898366195809, 0.824841025625018, 0.877962211974623, 0.9227853738926778, 0.9583291958881769, 0.9838388820037522], [0.001151677857036659, 0.007007259641384841, 0.0203148927103142, 0.0414025594490785, 0.06777171118707008, 0.09471990539100715, 0.11673118519090996, 0.1291364544780523, 0.12947566139904687, 0.11812235960758868, 0.09800772590075595, 0.07359961540231305, 0.04953526814944095, 0.02938957362654602, 0.0149639646714071, 0.00625271092938539, 0.001979666990405483, 0.00040351473706229986, 3.429268084906047e-5]))]), AbstractMeasure[Beta01Measure(PolyChaos.var"#102#103"{Float64, Float64}(2.1, 3.2), (0.0, 1.0), false, 2.1, 3.2), Measure("my_meas", Main.w, (0.0, 1.0), false, Dict{Any, Any}())]), AbstractOrthoPoly{M, Quad{Float64, Vector{Float64}}} where M<:AbstractMeasure[Beta01OrthoPoly{Vector{Float64}, Beta01Measure, Quad{Float64, Vector{Float64}}}(4, [0.39622641509433965, 0.45308865339881105, 0.4732655766681396, 0.482729089351984, 0.4879233481934926, 0.49108064278342917, 0.49314292190864784, 0.4945640779897571, 0.4955849084142159, 0.4963428640512603, 0.49692106736331404, 0.4973721930243669, 0.49773093800397555, 0.4980209139779084, 0.4982586420286101, 0.4984559630451982, 0.4986215434156863, 0.4987618443402984, 0.4988817625639983, 0.4989850639437675], [1.0, 0.03797318144060756, 0.049528193464418745, 0.054496397032519545, 0.05707494398498378, 0.05858224318366976, 0.05953877522858778, 0.06018346810103027, 0.0606384719605267, 0.06097150420824787, 0.061222558056958205, 0.06141648535528466, 0.061569388295544744, 0.061692071403646374, 0.06179200321800959, 0.061874480034347076, 0.06194334217794196, 0.06200142915302492, 0.06205087705273863, 0.062093317746215744], Beta01Measure(PolyChaos.var"#102#103"{Float64, Float64}(2.1, 3.2), (0.0, 1.0), false, 2.1, 3.2), Quad{Float64, Vector{Float64}}("golubwelsch", 19, [0.00878890525374891, 0.028450106116514506, 0.05852125536671723, 0.09833813997192814, 0.14702350691009355, 0.20350504892083368, 0.266538831062297, 0.3347366435399401, 0.40659656306145764, 0.48053602819819125, 0.5549266957500689, 0.6281303108338084, 0.6985348047519305, 0.7645898365922597, 0.8248410256009169, 0.8779622119533655, 0.9227853738738273, 0.9583291958712525, 0.9838388819882817], [0.00115167784765143, 0.007007259624729887, 0.020314892689531894, 0.04140255942833847, 0.06777171117046146, 0.09471990538153915, 0.11673118518989155, 0.12913645448496422, 0.12947566141190417, 0.11812235962361355, 0.09800772591714006, 0.07359961541684255, 0.04953526816082873, 0.029389573634464723, 0.014963964676256315, 0.006252710931950127, 0.0019796669915308576, 0.00040351473743713024, 3.429268092282879e-5])), OrthoPoly{Vector{Float64}, Measure, Quad{Float64, Vector{Float64}}}("my_op", 4, [0.39622641509448137, 0.45308865339922233, 0.47326557666900376, 0.4827290893535496, 0.48792334819610145, 0.4910806427875327, 0.4931429219148361, 0.4945640779987818, 0.495584908427019, 0.4963428640690031, 0.49692106738740605, 0.4973721930564992, 0.4977309380461525, 0.4980209140324845, 0.49825864209832377, 0.49845596313321483, 0.4986215435256319, 0.4987618444763114, 0.4988817627307632, 0.4989850641465735], [0.9999999999996537, 0.037973181440564815, 0.049528193464296905, 0.054496397032274234, 0.05707494398455889, 0.05858224318299423, 0.059538775227572845, 0.06018346809956606, 0.060638471958479156, 0.060971504205455586, 0.0612225580532287, 0.061416485350391216, 0.06156938828922262, 0.06169207139558928, 0.06179200320786561, 0.06187448002171521, 0.06194334216236786, 0.062001429133996595, 0.06205087702968264, 0.06209331771849251], Measure("my_meas", Main.w, (0.0, 1.0), false, Dict{Any, Any}()), Quad{Float64, Vector{Float64}}("golubwelsch", 19, [0.008788905307716177, 0.02845010617907043, 0.0585212554299101, 0.09833814003385198, 0.14702350696982944, 0.20350504897777927, 0.2665388311160098, 0.33473664359009403, 0.40659656310782666, 0.4805360282406439, 0.5549266957885687, 0.6281303108684106, 0.6985348047827785, 0.7645898366195809, 0.824841025625018, 0.877962211974623, 0.9227853738926778, 0.9583291958881769, 0.9838388820037522], [0.001151677857036659, 0.007007259641384841, 0.0203148927103142, 0.0414025594490785, 0.06777171118707008, 0.09471990539100715, 0.11673118519090996, 0.1291364544780523, 0.12947566139904687, 0.11812235960758868, 0.09800772590075595, 0.07359961540231305, 0.04953526814944095, 0.02938957362654602, 0.0149639646714071, 0.00625271092938539, 0.001979666990405483, 0.00040351473706229986, 3.429268084906047e-5]))])
mt2 = Tensor(2,mop)
mt3 = Tensor(3,mop)
mT2 = [ mt2.get([i,i]) for i=0:dim(mop)-1 ]
15-element Vector{Float64}:
 0.9999999999996537
 0.03797318144059441
 0.037973181440551666
 0.0018807430768492354
 0.00144196250871918
 0.0018807430768424916
 0.00010249372143213834
 7.141779810028215e-5
 7.141779810010645e-5
 0.00010249372143130948
 5.8498234095518276e-6
 3.892012680461295e-6
 3.5371945211048704e-6
 3.892012680434202e-6
 5.849823409460971e-6

Notice that mT2 carries the elements of the 2-dimensional tensors for the univariate bases opq and my_opq. The encoding is given by the multi-index mop.ind

mop.ind
15×2 Matrix{Int64}:
 0  0
 1  0
 0  1
 2  0
 1  1
 0  2
 3  0
 2  1
 1  2
 0  3
 4  0
 3  1
 2  2
 1  3
 0  4

To cross-check the results we can distribute the multi-index back to its univariate indices with the help of findUnivariateIndices.

ind_opq = findUnivariateIndices(1,mop.ind)
ind_my_opq = findUnivariateIndices(2,mop.ind)
5-element Vector{Int64}:
  1
  3
  6
 10
 15
mT2[ind_opq] - normsq
5-element Vector{Float64}:
 -3.462785613805863e-13
 -1.3149203947904198e-14
 -6.51171824794794e-16
 -3.54940686217442e-17
 -2.025255776885032e-18
mT2[ind_my_opq] - my_normsq
5-element Vector{Float64}:
 0.0
 0.0
 0.0
 0.0
 0.0