WIP univariate diffuse smoother

DynareJulia · Feb 13, 2023 · 934f8ef · 934f8ef
1 parent 2c2d4ee
commit 934f8ef
Show file tree

Hide file tree

Showing 3 changed files with 13 additions and 5 deletions.
diff --git a/doc/notes.tex b/doc/notes.tex
@@ -236,8 +236,7 @@ \section{Univariate diffuse smoother step}
 \end{align*}
 
 For $i = p_{t-1},\ldots,0$, if $|F_{t,i}| > 0$,
-\begin{align*}
-  
+\begin{align*}  
   r0_{t-1,i-1} &= L'_{\infty,i}r0_{t,i} \\
   r1_{t-1,i-1} &= Z'_{t,i}F^{-1}_{t,i}\nu_{t,i} +
                  L'_{\infty,t,i}r0_{t,i} + L'_{0,t,i}r1_{t,i} \\

diff --git a/src/univariate_step.jl b/src/univariate_step.jl
@@ -223,6 +223,7 @@ function diffuse_univariate_step!(Y, t, Z, H, T, RQR, a, Pinf, Pstar, diffuse_ka
         if Finf > diffuse_kalman_tol                 # F_{\infty,t,i} = 0, use upper part of bracket on p. 175 DK (2012) for w_{t,i}
             copy!(K0, ZPinf)
             rmul!(K0, 1/Finf)
+            @show K0
             a .+= v.*K0
             # Pstar     = Pstar + K0*(K0_Finf'*(Fstar/Finf)) - K1*K0_Finf' - K0_Finf*K1'
             ger!( Fstar/Finf, ZPinf, K0, Pstar)
@@ -401,6 +402,7 @@ function univariate_diffuse_smoother_step!(T::AbstractMatrix{W},
                                            ws::DiffuseKalmanSmootherWs) where {W <: AbstractFloat,  U <: Real}
     ny = size(Finf, 1)
     ns = size(L0, 1)
+    @show v
     @show r0_1
     @show r1_1
     for i = ny: -1: 1
@@ -418,12 +420,17 @@ function univariate_diffuse_smoother_step!(T::AbstractMatrix{W},
             # get_L1!(L1, K0, K1, Finf, Fstar, Z, i)  
             fill!(L1, 0.0)
             ger!(-1.0, vK1, vZ, L1)
+            @show vK1
+            @show vZ
             # compute r1_{t,i-1} first because it depends
             # upon r0{t,i} 
             # update_r1!(r1, r0, Z, v, Finv, L0, L1, i)
             copy!(ws.tmp_ns, vZ)
             rmul!(ws.tmp_ns, v[i]*iFinf)
+            @show ws.tmp_ns
+            @show L1
             mul!(ws.tmp_ns, transpose(L1), r0_1, 1.0, 1.0)
+            @show ws.tmp_ns
             mul!(ws.tmp_ns, transpose(L0), r1_1, 1.0, 1.0)
             copy!(r1, ws.tmp_ns)
             #  r0(:,t) = Linf'*r0(:,t);   % KD (2000), eq. (25) for r_0

diff --git a/test/test_univariate_step.jl b/test/test_univariate_step.jl
@@ -380,6 +380,8 @@ pstar0 = copy(Pstar[:, :, 1])
 pstar = copy(pstar0)
 QQ = R*Q*R'
 KalmanFilterTools.diffuse_univariate_step!(y, t, Z, H, T, QQ, a, pinf, pstar, 1e-10, 1e-10, ws_d)
+@show ws_d.K0[:, :, 1]
+@show ws_d.K[:, :, 1]
 v = y - c - Z*a0
 K0 = T*pinf0*Z'*inv(Z*pinf0*Z')
 a1 = T*a0 + K0*v
@@ -396,8 +398,8 @@ K = copy(K0)
 
 L0_target = I(ns) - K0*Z
 L1_target = -K1*Z
-r1_target = Z'inv(Finf)*v + L0_target'*r1 +L1_target'*r0
-r0_target = L0_target'*r0
+r1_target = Z'inv(Finf)*v + L0_target'*r1_1 +L1_target'*r0_1
+r0_target = L0_target'*r0_1
 N0_target = L0_target'*N0*L0_target
 N1_target = Z'*inv(Finf)*Z + L0_target'*N1*L0_target + L1_target'*N0_target*L0_target
 Fstar = -inv(Finf)*Fstar*inv(Finf)
@@ -407,7 +409,7 @@ tol = 1e-12
 KalmanFilterTools.univariate_diffuse_smoother_step!(T, ws_d.F[:, :, 1], ws_d.Fstar[:, :, 1],
 ws_d.K0[:, :, 1], ws_d.K[:, :, 1],
 ws_d.L, ws_d.L1, ws_d.N, ws_d.N1,
-ws_d.N2, r0, r1, ws_d.v[:,1], Z,
+ws_d.N2, r0, r1, r0_1, r1_1, ws_d.v[:,1], Z,
 pinf, pstar, tol, ws_d)
 
 #=