Kicking the energy at every turn using the result of the normal form

The normal form normal_spin 6.2 contains a flag nf%stochastic. If true, FPP computes a matrix nf%STOCH(6,6) and 3 stochastic kicks in nf%KICK(3). The map is gotten by the following formula:

$$\mbox{\Huge${\vec{x}}_{n+1} = $} \mbox{\Huge${S}^{-1}$}\left({\mbox{\Huge$SM{\vec{x}}_{n}$}+\left(... ...}_{2}\cr {t}_{2}{k}_{2}\cr {t}_{3}{k}_{3}\cr {t}_{3}{k}_{3}\cr}}\right)}\right)$$

$${\rm w}{\rm h}{\rm e}{\rm r}{\rm e} S= nf\%STOCH {\rm a}{\rm n}{\rm d} k= nf\%KICK$$ (6.4)

do i=1,n
 x=matmul(mat,x)
 x=matmul(stoch,x)
 
do j=1,3
       t=RANF()
       if(t>half) then
          t=one
       else
          t=-one
       endif
   x(2*j-1)= x(2*j-1)+ t*nf%kick(j)
   x(2*j   )= x(2*j)+  t*nf%kick(j)
 enddo
   x=matmul(stochi,x)
 
  sizes(1,1)=sizes(1,1)+x(1)**2
  sizes(2,2)=sizes(2,2)+x(2)**2
  sizes(1,3)=sizes(1,3)+x(1)*x(3)
enddo

sizes=sizes/n
write(6,*)  sizes(1,1), sizes(2,2), sizes(1,3)
write(6,*)  nf%s_ij0(1,1), nf%s_ij0(2,2), nf%s_ij0(1,3)