agd.Metrics.asym_rander

 1# Copyright 2020 Jean-Marie Mirebeau, University Paris-Sud, CNRS, University Paris-Saclay
 2# Distributed WITHOUT ANY WARRANTY. Licensed under the Apache License, Version 2.0, see http://www.apache.org/licenses/LICENSE-2.0
 3
 4import numpy as np
 5from . import misc
 6from .base import Base
 7from .riemann import Riemann
 8from .rander import Rander
 9from .asym_quad import AsymQuad
10
11from .. import LinearParallel as lp
12from .. import AutomaticDifferentiation as ad
13from .. import FiniteDifferences as fd
14
15class AsymRander(Base):
16	r"""
17	Asymmetric Rander norms take the form
18	$$
19	F(x) = \sqrt{<x,m x> + max(0,<u,x>)^2 + max(0,<v,x>)^2} + <w,x>
20	$$
21	where m is a given symmetric positive definite tensor, and 
22	$u,v,w$ are given vectors. The vector $w$ must be small enough, so that $F(x)>0$ 
23	for all $x\neq 0$.
24
25	Asymmetric Rander norms generalize both Rander norms and Asymmetric quadratic norms.
26	They were proposed by Da Chen in the context of image processing applications.
27
28	Member fields and __init__ arguments : 
29	- m : an array of shape (vdim,vdim,n1,..,nk) where vdim is the ambient space dimension.
30	The array must be symmetric, a.k.a m[i,j] = m[j,i] for all 0<=i<j<vdim.
31	- u,v,w : an array of shape (vdim,n1,...,nk)
32	"""
33
34	def __init__(self,m,u,v,w):
35		m = ad.asarray(m)
36		if u is None: u=np.zeros_like(m[0])
37		if v is None: v=np.zeros_like(m[0])
38		if w is None: w=np.zeros_like(m[0])
39		m,u,v,w = (ad.asarray(e) for e in (m,u,v,w))
40		self.m,self.u,self.v,self.w = fd.common_field((m,u,v,w),(2,1,1,1))
41
42	def norm(self,x):
43		x,m,u,v,w = fd.common_field((ad.asarray(x),self.m,self.u,self.v,self.w),(1,2,1,1,1))
44		return np.sqrt(lp.dot_VAV(x,m,x) + np.maximum(lp.dot_VV(x,u),0.)**2
45			+ np.maximum(lp.dot_VV(x,v),0.)**2 ) + lp.dot_VV(x,w)
46
47	@property
48	def vdim(self): return len(self.m)
49
50	@property
51	def shape(self): return self.m.shape[2:]	
52
53	def flatten(self):
54		m,u,v,w = self.m,self.u,self.v,self.w
55		return np.concatenate((misc.flatten_symmetric_matrix(m),u,v,w),axis=0)
56
57	def model_HFM(self):
58		return "AsymRander"+str(self.vdim)
59
60	def __iter__(self):
61		yield self.m
62		yield self.u
63		yield self.v
64		yield self.w
65
66	@classmethod
67	def from_cast(cls,metric):
68		if isinstance(metric,cls): return metric
69		elif isinstance(metric,Rander): 
70			z = np.zeros_like(metric.w)
71			return cls(metric.m,z,z,metric.w)
72		else:
73			asym = AsymQuad.from_cast(metric)
74			z = np.zeros_like(metric.w)
75			return cls(metric.m,metric.w,z,z)
76
77	def make_proj_dual(self,**kwargs):
78		"""
79		Orthogonal projection onto the dual unit ball.
80		- **kwargs : passed to AsymQuad.make_proj_dual
81		"""
82		if not np.allclose(v,0.): raise ValueError("Sorry, prox compuation requires v=0")
83		proj_mu = AsymQuad(self.m,self.u).make_proj_dual(**kwargs)
84		return lambda x : proj_mu(x-self.w)+self.w
class AsymRander(agd.Metrics.base.Base):
16class AsymRander(Base):
17	r"""
18	Asymmetric Rander norms take the form
19	$$
20	F(x) = \sqrt{<x,m x> + max(0,<u,x>)^2 + max(0,<v,x>)^2} + <w,x>
21	$$
22	where m is a given symmetric positive definite tensor, and 
23	$u,v,w$ are given vectors. The vector $w$ must be small enough, so that $F(x)>0$ 
24	for all $x\neq 0$.
25
26	Asymmetric Rander norms generalize both Rander norms and Asymmetric quadratic norms.
27	They were proposed by Da Chen in the context of image processing applications.
28
29	Member fields and __init__ arguments : 
30	- m : an array of shape (vdim,vdim,n1,..,nk) where vdim is the ambient space dimension.
31	The array must be symmetric, a.k.a m[i,j] = m[j,i] for all 0<=i<j<vdim.
32	- u,v,w : an array of shape (vdim,n1,...,nk)
33	"""
34
35	def __init__(self,m,u,v,w):
36		m = ad.asarray(m)
37		if u is None: u=np.zeros_like(m[0])
38		if v is None: v=np.zeros_like(m[0])
39		if w is None: w=np.zeros_like(m[0])
40		m,u,v,w = (ad.asarray(e) for e in (m,u,v,w))
41		self.m,self.u,self.v,self.w = fd.common_field((m,u,v,w),(2,1,1,1))
42
43	def norm(self,x):
44		x,m,u,v,w = fd.common_field((ad.asarray(x),self.m,self.u,self.v,self.w),(1,2,1,1,1))
45		return np.sqrt(lp.dot_VAV(x,m,x) + np.maximum(lp.dot_VV(x,u),0.)**2
46			+ np.maximum(lp.dot_VV(x,v),0.)**2 ) + lp.dot_VV(x,w)
47
48	@property
49	def vdim(self): return len(self.m)
50
51	@property
52	def shape(self): return self.m.shape[2:]	
53
54	def flatten(self):
55		m,u,v,w = self.m,self.u,self.v,self.w
56		return np.concatenate((misc.flatten_symmetric_matrix(m),u,v,w),axis=0)
57
58	def model_HFM(self):
59		return "AsymRander"+str(self.vdim)
60
61	def __iter__(self):
62		yield self.m
63		yield self.u
64		yield self.v
65		yield self.w
66
67	@classmethod
68	def from_cast(cls,metric):
69		if isinstance(metric,cls): return metric
70		elif isinstance(metric,Rander): 
71			z = np.zeros_like(metric.w)
72			return cls(metric.m,z,z,metric.w)
73		else:
74			asym = AsymQuad.from_cast(metric)
75			z = np.zeros_like(metric.w)
76			return cls(metric.m,metric.w,z,z)
77
78	def make_proj_dual(self,**kwargs):
79		"""
80		Orthogonal projection onto the dual unit ball.
81		- **kwargs : passed to AsymQuad.make_proj_dual
82		"""
83		if not np.allclose(v,0.): raise ValueError("Sorry, prox compuation requires v=0")
84		proj_mu = AsymQuad(self.m,self.u).make_proj_dual(**kwargs)
85		return lambda x : proj_mu(x-self.w)+self.w

Asymmetric Rander norms take the form $$ F(x) = \sqrt{ + max(0,)^2 + max(0,)^2} + $$ where m is a given symmetric positive definite tensor, and $u,v,w$ are given vectors. The vector $w$ must be small enough, so that $F(x)>0$ for all $x\neq 0$.

Asymmetric Rander norms generalize both Rander norms and Asymmetric quadratic norms. They were proposed by Da Chen in the context of image processing applications.

Member fields and __init__ arguments :

  • m : an array of shape (vdim,vdim,n1,..,nk) where vdim is the ambient space dimension. The array must be symmetric, a.k.a m[i,j] = m[j,i] for all 0<=i
  • u,v,w : an array of shape (vdim,n1,...,nk)
AsymRander(m, u, v, w)
35	def __init__(self,m,u,v,w):
36		m = ad.asarray(m)
37		if u is None: u=np.zeros_like(m[0])
38		if v is None: v=np.zeros_like(m[0])
39		if w is None: w=np.zeros_like(m[0])
40		m,u,v,w = (ad.asarray(e) for e in (m,u,v,w))
41		self.m,self.u,self.v,self.w = fd.common_field((m,u,v,w),(2,1,1,1))
def norm(self, x):
43	def norm(self,x):
44		x,m,u,v,w = fd.common_field((ad.asarray(x),self.m,self.u,self.v,self.w),(1,2,1,1,1))
45		return np.sqrt(lp.dot_VAV(x,m,x) + np.maximum(lp.dot_VV(x,u),0.)**2
46			+ np.maximum(lp.dot_VV(x,v),0.)**2 ) + lp.dot_VV(x,w)

Norm or quasi-norm defined by the class, often denoted $N$ in mathematical formulas. Unless incorrect data is provided, this member function obeys, for all vectors $u,v\in R^d$ and all $\alpha \geq 0$

  • $N(u+v) \leq N(u)+N(v)$
  • $N(\alpha u) = \alpha N(u)$
  • $N(u)\geq 0$ with equality iff $u=0$.

Broadcasting will occur depending on the shape of $v$ and of the class data.

vdim
48	@property
49	def vdim(self): return len(self.m)

The ambient vector space dimension, often denoted $d$ in mathematical formulas.

shape
51	@property
52	def shape(self): return self.m.shape[2:]	

Dimensions of the underlying domain. Expected to be the empty tuple, or a tuple of length vdim.

def flatten(self):
54	def flatten(self):
55		m,u,v,w = self.m,self.u,self.v,self.w
56		return np.concatenate((misc.flatten_symmetric_matrix(m),u,v,w),axis=0)

Flattens and concatenate the member fields into a single array.

def model_HFM(self):
58	def model_HFM(self):
59		return "AsymRander"+str(self.vdim)

The name of the 'model' for parameter, as input to the HFM library.

@classmethod
def from_cast(cls, metric):
67	@classmethod
68	def from_cast(cls,metric):
69		if isinstance(metric,cls): return metric
70		elif isinstance(metric,Rander): 
71			z = np.zeros_like(metric.w)
72			return cls(metric.m,z,z,metric.w)
73		else:
74			asym = AsymQuad.from_cast(metric)
75			z = np.zeros_like(metric.w)
76			return cls(metric.m,metric.w,z,z)

Produces a metric by casting another metric of a compatible type.

def make_proj_dual(self, **kwargs):
78	def make_proj_dual(self,**kwargs):
79		"""
80		Orthogonal projection onto the dual unit ball.
81		- **kwargs : passed to AsymQuad.make_proj_dual
82		"""
83		if not np.allclose(v,0.): raise ValueError("Sorry, prox compuation requires v=0")
84		proj_mu = AsymQuad(self.m,self.u).make_proj_dual(**kwargs)
85		return lambda x : proj_mu(x-self.w)+self.w

Orthogonal projection onto the dual unit ball.

  • **kwargs : passed to AsymQuad.make_proj_dual