Commit 5dc72793 authored by alexander.dornheim's avatar alexander.dornheim
Browse files

Started proting QNAP, not yet complete

parent 2613322c
# QNAP.m
# written on May 23, 2017 by
# Russell Luke
# Inst. Fuer Numerische und Angewandte Mathematik
# Universitaet Gottingen
#
# DESCRIPTION: Averaged prox mappings with a Quasi-Newton acceleration
#
# INPUT: method_input, a data structure
#
# OUTPUT: method_output, a data structure with
# u = the algorithm fixed point
# stats = [iter, change, gap] where
# gap = squared gap distance normalized by
# the magnitude constraint.
# change = the percentage change in the norm
# squared difference in the iterates
# iter = the number of iterations the algorithm
# performed
#
# USAGE: method_output = QNAP(method_input)
#
# Nonstandard Matlab function calls: samsara, self.prox1 and .Prox2
# It is assumed that the Prox1 is a projection onto the diagonal of the
# product space, P_Diag. Must have samsara, a reverse communication nonlinear
# optimization package installed
#
#
class QNAP(Algortihm):
def __init__(self,config):
self.samsare = samsara()
def run(self, u, tol, maxiter):
while iter < maxiter and change[iter] >= tol:
iter += 1
config['iter'] =iter
tmp_u = prox1.work(tmp1)
# make call to SAMSARA for acceleration step
# since for the product space Prox1 is the
# projection onto the diagonal, all arrays in
# tmp_u are the same, so we can just send one
# array to samsara. Have to reshape the complex
# array as a real-valued vector. Use the gap as the
# objective function, and the negative gradient tmp_u- u
if dim_number == 2:
if (isreal(u)):
if(iter>4):
uold_vec=u[:,1]
unew_vec=tmp_u[:,1]
gradfold_vec = gradfnew_vec
gradfnew_vec = (u[:,1]- tmp_u[:,1)]
unew_vec, uold_vec, gap[iter-1], gradfold_vec, change[iter], samsara_data =...
self.samsara(uold_vec, unew_vec, ...
gap[iter-2]*self.Nx*self.Ny,...
gap[iter-1]*self.Nx*self.Ny,...
gradfold_vec, gradfnew_vec, samsara_data)
for in range self.product_space_dimension:
tmp_u[:,j]=unew_vec
gap[iter-1]=gap[iter-1]/(self.Nx*self.Ny)
else:
unew_vec=tmp_u[:,1]
uold_vec=u[:,1]
gradfnew_vec = (u[:,1]- tmp_u[:,1])
else:
if iter>3:
uold_vec=reshape([real(u(:,1)), imag(u(:,1))], ...
self.Ny*self.Nx*2, 1);
unew_vec=reshape([real(tmp_u(:,1)), imag(tmp_u(:,1))], ...
self.Ny*self.Nx*2, 1)
gradfold_vec = gradfnew_vec
gradfnew_vec = uold_vec-unew_vec
unew_vec,uold_vec,gap[iter-1],gradfold_vec, change[iter), samsara_data=...
feval('samsara', uold_vec, unew_vec, ...
gap(iter-2)*self.Nx*self.Ny,...
gap(iter-1)*self.Nx*self.Ny,...
gradfold_vec, gradfnew_vec, samsara_data);
# now reshape unew_vec
tmp_u_vec = unew_vec(1:self.Ny*self.Nx)+1i*unew_vec(self.Ny*self.Nx+1:self.Ny*self.Nx*2);
for j in range(self.product_space_dimension)
tmp_u[:,j]=tmp_u_vec;
gap[iter-1]=gap[iter-1]/(self.Nx*self.Ny)
else:
uold_vec=reshape([real(u(:,1)), imag(u(:,1))], ...
self.Ny*self.Nx*2, 1);
unew_vec=reshape([real(tmp_u(:,1)), imag(tmp_u(:,1))], ...
self.Ny*self.Nx*2, 1);
gradfnew_vec = uold_vec-unew_vec
elif dim_number == 3
if (isreal(u))
tmp2 = u[:,:,1]
uold_vec=reshape(tmp2,self.Nx*self.Ny,1)
tmp2 = tmp_u[:,:,1]
unew_vec = reshape(tmp2,self.Nx*self.Ny,1);
if iter<=3:
gradfnew_vec = uold_vec-unew_vec
else:
gradfold_vec = gradfnew_vec
gradfnew_vec = uold_vec-unew_vec
unew_vec,uold_vec,gap[iter-1],gradfold_vec, change[iter), samsara_data=...
feval('samsara', uold_vec, unew_vec, ...
gap(iter-2)*self.Nx*self.Ny,...
gap(iter-1)*self.Nx*self.Ny,...
gradfold_vec, gradfnew_vec, samsara_data);
# now reshape and replace u
gap[iter-1]=gap[iter-1]/(self.Nx*self.Ny);
tmp2=reshape(unew_vec,self.Ny,self.Nx);
for j in range(self.product_space_dimension):
tmp_u[:,:,j]=tmp2
end
else:
tmp2=[real(u[:,:,1]),imag(u[:,:,1)]]
uold_vec=reshape(tmp2,self.Nx*self.Ny*2,1)
tmp2=[real(tmp_u(:,:,1)), imag(tmp_u(:,:,1))]
unew_vec = reshape(tmp2,self.Nx*self.Ny*2,1)
if iter<=3:
gradfnew_vec = uold_vec-unew_vec
else
gradfold_vec = gradfnew_vec
gradfnew_vec = uold_vec-unew_vec
unew_vec,uold_vec,gap[iter-1],gradfold_vec, change[iter), samsara_data=...
feval('samsara', uold_vec, unew_vec, ...
gap(iter-2)*self.Nx*self.Ny,...
gap(iter-1)*self.Nx*self.Ny,...
gradfold_vec, gradfnew_vec, samsara_data);
# now reshape and replace u
gap[iter-1]=gap[iter-1]/(self.Nx*self.Ny);
tmp2=reshape(unew_vec,self.Ny,self.Nx*2);
unew=tmp2(:,1:self.Nx) +1i*tmp2(:,self.Nx+1:2*self.Nx);
for j=1:self.product_space_dimension
tmp_u[:,:,j]=unew
else: # product space of 3D arrays
print('Cannot handle 3D arrays on the product space yet')
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment