From e951e0c7eb1553bcbaf9a00bc4565b7bd5424a5d Mon Sep 17 00:00:00 2001
From: Yuqing <yuqingzhu@ucsb.edu>
Date: Thu, 21 Apr 2022 01:21:52 -0700
Subject: [PATCH 1/8] Update implementations of svt paper

---
 autodp/dp_bank.py                      | 17 +++++++++--
 autodp/mechanism_zoo.py                |  8 +++---
 autodp/rdp_bank.py                     | 40 +++++++++++++++-----------
 papers/NeurIPS20-SVT/exp2_align_var.py | 19 ++++++------
 4 files changed, 53 insertions(+), 31 deletions(-)

diff --git a/autodp/dp_bank.py b/autodp/dp_bank.py
index 4fbcbaf..bcc0028 100644
--- a/autodp/dp_bank.py
+++ b/autodp/dp_bank.py
@@ -85,22 +85,33 @@ def fun(x):
     else:
         raise RuntimeError(f"Failed to find epsilon: {results.flag}")
 
-def eps_generalized_gaussian(x, sigma, delta,k, c, c_tilde):
+def eps_generalized_gaussian(params, delta):
     """
     submodule for generalized SVT with Gaussian noise
     we want to partition c into [c/c'] parts, each part using (k choose c')
     need to check whether (k choose c') > log(1/delta')
     k is the maximam number of queries to answer for each chunk
-    x is log delta for each chunk, it needs to be negative
     :param x:
     :param sigma:
     :param delta:
     :return:
     """
+    sigma = params['sigma']
+    k = params['k']
+    c = params['c']
+    # delta' (delta1) is the delta budget for each chunk.
+    c_tilde = int(np.sqrt(c))
+    # divide half of the delta budget into c_tilde chunks, x is log delta for each chunk, it needs to be negative.
+    x = np.log(delta/(2 * c_tilde))
     acct = dp_acct.DP_acct()
-    per_delta = np.exp(x) # per_delta for each c' chunk
     coeff = comb(k,c_tilde)
+    # Find a proper per delta such that c_tilde * per_delta < delta and (k choose c') > log(1/per_delta)
+    per_delta = min( 0.5 /coeff, delta/(2*c_tilde))
+    x = np.log(per_delta)
+    if per_delta > 1.0/(coeff):
+        print('per delta', per_delta, 'coeff', coeff, 'k', k, 'c', c)
     assert per_delta < 1.0/(coeff)
+
     #compute the eps per step with 1/(sigma_1**2) + sqrt(2/simga_1**2 *(log k + log(1/epr_delta)))
     # compose eps for each chunk
     while c:
diff --git a/autodp/mechanism_zoo.py b/autodp/mechanism_zoo.py
index e9d9532..e10ef44 100644
--- a/autodp/mechanism_zoo.py
+++ b/autodp/mechanism_zoo.py
@@ -347,11 +347,11 @@ def __init__(self,params,name='GaussianSVT', rdp_c_1=True):
         self.name=name
         if rdp_c_1 == True:
             self.name = name + 'c_1'
-            self.params = {'sigma': params['sigma'], 'k': params['k'], 'margin':params['margin']}
+            self.params = {'sigma': params['sigma'], 'sigma_nu': params['sigma_nu'], 'k': params['k'], 'margin':params['margin']}
             new_rdp = lambda x: rdp_bank.RDP_gaussian_svt_c1(self.params, x)
         else:
             self.name = name + 'c>1'
-            self.params = {'sigma':params['sigma'],'k':params['k'], 'c':params['c']}
+            self.params = {'sigma':params['sigma'],'sigma_nu': params['sigma_nu'], 'k':params['k'], 'c':params['c']}
             new_rdp = lambda x: rdp_bank.RDP_gaussian_svt_cgreater1(self.params, x)
         self.propagate_updates(new_rdp, 'RDP')
 
@@ -384,11 +384,11 @@ def __init__(self, params=None,approxDP_off=False, name='StageWiseMechanism'):
         Mechanism.__init__(self)
 
         self.name = name # When composing
-        self.params = {'sigma': params['sigma'], 'k':params['k'], 'c':params['c']}
+        self.params = params
         self.delta0 = 0
 
         if not approxDP_off:  # Direct implementation of approxDP
-            new_approxdp = lambda x: dp_bank.get_generalized_gaussian(params, x)
+            new_approxdp = lambda x: dp_bank.eps_generalized_gaussian(params, x)
             self.propagate_updates(new_approxdp, 'approxDP_func')
 
 
diff --git a/autodp/rdp_bank.py b/autodp/rdp_bank.py
index f7eaacc..e0b46ca 100644
--- a/autodp/rdp_bank.py
+++ b/autodp/rdp_bank.py
@@ -12,7 +12,7 @@
 import numpy as np
 import math
 from autodp import utils
-from autodp.utils import _log1mexp
+from autodp.utils import _log1mexp, logcomb, stable_logsumexp
 
 
 def stable_log_diff_exp(x):
@@ -359,15 +359,20 @@ def RDP_gaussian_svt_cgreater1(params, alpha):
 
     Args:
         k:the maximum length before svt stops
-        sigma: noise added to the threshold.
+        sigma: noise added to the threshold divide by L2 sensitivity. 
+		sigma_nu: noise added to the query divide by 2 * L2 sensitivity.
         c: the cut-off parameter in SVT.
     """
-    sigma = params['sigma']
+    sigma_rho = params['sigma']
+    sigma_nu = params['sigma_nu']
     c = max(params['c'], 1)
     k = params['k']  # the algorithm stops either k is achieved or c is achieved
-    rdp_rho = 0.5 / (sigma ** 2) * alpha
-    c_log_n_c = c * np.log(k / c)
-    ret_rdp = c_log_n_c * 1.0 / (alpha - 1) + rdp_rho * (c + 1)
+    rdp_rho = 0.5 / (sigma_rho ** 2) * alpha
+    rdp_nu = 0.5 / (sigma_nu ** 2) * alpha
+    #c_log_n_c = c * np.log(k / c) # approximate bound
+    log_comb = [logcomb(k, i) for i in range(c+1)]
+    log_sum_comb = stable_logsumexp(log_comb)
+    ret_rdp = log_sum_comb * 1.0 / (alpha - 1) + rdp_rho + rdp_nu * c
     return ret_rdp
 
 
@@ -380,27 +385,30 @@ def RDP_gaussian_svt_c1(params, alpha):
 
     Args:
         k:the maximum length before svt stops
-        sigma: noise added to the threshold.
+        sigma: noise added to the threshold, divided by L2 sensitivity.
+		sigma_nu: noise added to the query, divided by the 2*L2 sensitivity.
         c: the cut-off parameter in SVT.
     """
     sigma = params['sigma']
+    sigma_nu = params['sigma_nu']
     k = params['k']
     margin = params['margin']
-    c = 1
-
-
+    
     rdp_rho = 0.5 / (sigma ** 2) * alpha
-
-    ret_rdp = np.log(k) / (alpha - 1) + rdp_rho * 2
+    rdp_nu = 0.5 / (sigma **2 ) * alpha
+    ret_rdp = np.log(k) / (alpha - 1) + rdp_rho + rdp_nu
     if alpha == 1:
-        return ret_rdp * c
-    ################ Implement corollary 15 in NeurIPS-20
+        return ret_rdp
+    if sigma * np.sqrt(3) > sigma_nu:
+        return ret_rdp
+    # The code below Implements Proposition 10 in NeurIPS-20 with gamma=2, which requires
+	# (1) sigma_nu> sqrt{3}sigma_rho (2) queries are non-negative.
     inside_part = np.log(2 * np.sqrt(3) * math.pi * (1 + 9 * margin ** 2 / (sigma ** 2)))
     moment_term = utils.stable_logsumexp_two(0, inside_part + margin ** 2 * 1.0 / (sigma ** 2))
     moment_term = moment_term / (2.0 * (alpha - 1))
-    moment_based = moment_term + rdp_rho * 2
+    moment_based = moment_term + rdp_rho * 2 + rdp_nu 
 
-    return min(moment_based, ret_rdp) * c
+    return min(moment_based, ret_rdp) 
 
 
 
diff --git a/papers/NeurIPS20-SVT/exp2_align_var.py b/papers/NeurIPS20-SVT/exp2_align_var.py
index cd17a3b..ddd6cb5 100644
--- a/papers/NeurIPS20-SVT/exp2_align_var.py
+++ b/papers/NeurIPS20-SVT/exp2_align_var.py
@@ -7,7 +7,7 @@
 import matplotlib.pyplot as plt
 from scipy.stats import norm, laplace
 from scipy.special import comb
-import matplotlib.font_manager as fm
+#import matplotlib.font_manager as fm
 from autodp.mechanism_zoo import ExactGaussianMechanism, PureDP_Mechanism,SubsampleGaussianMechanism, GaussianMechanism, ComposedGaussianMechanism,GaussianSVT_Mechanism, NoisyScreenMechanism
 from autodp.transformer_zoo import Composition, AmplificationBySampling
 
@@ -26,6 +26,7 @@
 delta = 1e-6
 lambda_rho = 120
 lambda_nu = 240
+# Under this setting, sigma_1 / L2 sensitivity = sigma_2 /2 * L2 sensitivity
 sigma_1 = lambda_rho*np.sqrt(2)
 sigma_2 = 2*sigma_1
 eps_1 = 1.0 / lambda_rho
@@ -102,11 +103,12 @@ def exp_2a():
 
             # stage-wise generalized SVT, k is the maximum length of each chunk
             k = int(idx / np.sqrt(count_gau))
-            generalized_mech = StageWiseMechanism({'sigma':sigma_1,'k':k, 'c':count_gau})
+            # assume sensitivity is 1, sigma_nu denotes noise add to query / 2
+            generalized_mech = StageWiseMechanism({'sigma':sigma_1,'k':k, 'c':count_gau, 'delta':delta})
             eps_kov.append(generalized_mech.get_approxDP(delta))
 
             # Gaussian-SVT c>1 with RDP, k is the total length before algorithm stops
-            gaussianSVT_c = GaussianSVT_Mechanism({'sigma':sigma_1,'k':idx, 'c':count_gau}, rdp_c_1=False)
+            gaussianSVT_c = GaussianSVT_Mechanism({'sigma':sigma_1,'sigma_nu':sigma_1,'k':idx, 'c':count_gau}, rdp_c_1=False)
             eps_g.append(gaussianSVT_c.get_approxDP(delta))
 
             #Gaussian-SVT with c=1, we use average_k as the approximate maximum length of each chunk, margin is used in Proposition 10
@@ -114,6 +116,7 @@ def exp_2a():
             params_SVT = {}
             params_SVT['k'] = average_k
             params_SVT['sigma'] = sigma_1
+            params_SVT['sigma_nu'] = sigma_1
             params_SVT['margin'] = margin
             per_gaussianSVT_mech = GaussianSVT_Mechanism(params_SVT)
             gaussianSVT_mech = compose([per_gaussianSVT_mech],[max(count_gau, 1)])
@@ -129,7 +132,7 @@ def exp_2a():
             'weight': 'bold',
             'size': 18}
 
-    props = fm.FontProperties(family='Gill Sans', fname='/Library/Fonts/GillSans.ttc')
+    #props = fm.FontProperties(family='Gill Sans', fname='/Library/Fonts/GillSans.ttc')
     f, ax = plt.subplots()
     plt.figure(num=0, figsize=(12, 8), dpi=80, facecolor='w', edgecolor='k')
     plt.loglog(index, eps_a, '-r', linewidth=2)
@@ -148,7 +151,7 @@ def exp_2a():
     plt.yticks(fontsize=20)
     plt.xlabel(r'Iterations', fontsize=20)
     plt.ylabel(r'$\epsilon$', fontsize=20)
-    ax.set_title('Title', fontproperties=props)
-
-    plt.savefig('exp2a.pdf', bbox_inches='tight')
-exp_2a()
\ No newline at end of file
+    #ax.set_title('Title', fontproperties=props)
+    plt.show()
+    #plt.savefig('exp2a.pdf', bbox_inches='tight')
+exp_2a()

From 5096bc9b3397b5bd59f69ab775c2ee13bca788b5 Mon Sep 17 00:00:00 2001
From: Jonathan Lebensold <lebensold@meta.com>
Date: Fri, 7 Jul 2023 15:16:50 -0400
Subject: [PATCH 2/8] Support for global sensitivity. Fixed division by zero
 for gaussian SVT

---
 autodp/mechanism_zoo.py |  3 +-
 autodp/rdp_bank.py      | 79 ++++++++++++++++++++++++++++-------------
 2 files changed, 56 insertions(+), 26 deletions(-)

diff --git a/autodp/mechanism_zoo.py b/autodp/mechanism_zoo.py
index e10ef44..12a00ec 100644
--- a/autodp/mechanism_zoo.py
+++ b/autodp/mechanism_zoo.py
@@ -347,7 +347,8 @@ def __init__(self,params,name='GaussianSVT', rdp_c_1=True):
         self.name=name
         if rdp_c_1 == True:
             self.name = name + 'c_1'
-            self.params = {'sigma': params['sigma'], 'sigma_nu': params['sigma_nu'], 'k': params['k'], 'margin':params['margin']}
+            valid_keys = ['sigma', 'k', 'margin', 'sigma_nu', 'Delta', 'sigma2', 'gamma']
+            self.params = dict(filter(lambda tuple: tuple[0] in valid_keys, params.items()))
             new_rdp = lambda x: rdp_bank.RDP_gaussian_svt_c1(self.params, x)
         else:
             self.name = name + 'c>1'
diff --git a/autodp/rdp_bank.py b/autodp/rdp_bank.py
index e0b46ca..f81e107 100644
--- a/autodp/rdp_bank.py
+++ b/autodp/rdp_bank.py
@@ -381,34 +381,63 @@ def RDP_gaussian_svt_c1(params, alpha):
     The RDP of the gaussian-based SVT with the cut-off parameter c=1.
 
     The detailed algorithm is described in Theorem 8 in
-    https://papers.nips.cc/paper/2020/file/e9bf14a419d77534105016f5ec122d62-Paper.pdf/
+    https://papers.nips.cc/paper/2020/file/e9bf14a419d77534105016f5ec122d62-Paper.pdf
 
-    Args:
-        k:the maximum length before svt stops
-        sigma: noise added to the threshold, divided by L2 sensitivity.
-		sigma_nu: noise added to the query, divided by the 2*L2 sensitivity.
-        c: the cut-off parameter in SVT.
-    """
-    sigma = params['sigma']
-    sigma_nu = params['sigma_nu']
-    k = params['k']
-    margin = params['margin']
+    A special case, where $k$ is not required, occurs when:
+    1. Queries are non-negative, and have low-sensitivity (assumed),
+    2. $\sigma_2 \geq \sqrt{3} \sigma_1$,
+    3. $\gamma = 2$
     
-    rdp_rho = 0.5 / (sigma ** 2) * alpha
-    rdp_nu = 0.5 / (sigma **2 ) * alpha
-    ret_rdp = np.log(k) / (alpha - 1) + rdp_rho + rdp_nu
-    if alpha == 1:
-        return ret_rdp
-    if sigma * np.sqrt(3) > sigma_nu:
+    Calculates in the standard bound and the special case (Proposition 10) when applicable.
+
+    params args:
+        k (int, optional): the maximum length before svt stops (used in Theorem 8. -- general case)
+        sigma (float): noise added to the threshold ($\sigma_1$ ($\sigma_{\\rho}$) in [ZW2020]).
+        sigma_nu (float): noise added to each query ($\sigma_2$ ($\sigma_{\nu}$) in [ZW2020]).
+            Defaults to \sqrt{3} \sigma_{1}
+        Delta (float): global sensitivity for each query
+        margin (float): query threshold parameter ($T$ in [ZW2020])
+        gamma (float): positive multiplicative factor (default=2)
+    """
+    # API expects `sigma` to be defined for all gaussian mechanisms, let sigma=sigma_1:
+    sigma_1 = params['sigma']
+    # [ZW2020] denotes sigma_nu as sigma_2. By default, calibrate query noise to apply Prop. 10:
+    sigma_2 = params.get('sigma_nu', np.sqrt(3) * sigma_1)
+    # AutoDP uses sigma2 as a param, for consistency we support it
+    sigma_2 = params.get('sigma2', sigma_2)
+    margin = params.get('margin', 0.)
+    Delta = params.get('Delta', 1.)
+    gamma = params.get('gamma', 2.)
+
+    rdp_rho = alpha * Delta**2 / (sigma_1 ** 2)
+    # By Theorem 8, $\eps_{\nu}(\alpha)$ is upper bounded by 2 \Delta.
+    rdp_nu = 2 * alpha * Delta**2 / (sigma_2 **2 )
+
+    # If $k$ provided, compute the standard bound
+    ret_rdp = np.inf
+    if 'k' in params:
+        k = params['k']
+        ret_rdp = np.divide(np.log(k) , alpha - 1, where=alpha!=1) + rdp_rho + rdp_nu
+    
+    # Check if the conditions for the Propisition 10 are satisfied:
+    if sigma_2 < np.sqrt(gamma + 1) * sigma_1:
+        # Check if k is defined, if not, no bound can be applied correctly
+        assert 'k' in params, "$k$ must be defined if Proposition 10 is not applicable."
         return ret_rdp
-    # The code below Implements Proposition 10 in NeurIPS-20 with gamma=2, which requires
-	# (1) sigma_nu> sqrt{3}sigma_rho (2) queries are non-negative.
-    inside_part = np.log(2 * np.sqrt(3) * math.pi * (1 + 9 * margin ** 2 / (sigma ** 2)))
-    moment_term = utils.stable_logsumexp_two(0, inside_part + margin ** 2 * 1.0 / (sigma ** 2))
-    moment_term = moment_term / (2.0 * (alpha - 1))
-    moment_based = moment_term + rdp_rho * 2 + rdp_nu 
-
-    return min(moment_based, ret_rdp) 
+
+    # Proposition 10. Bound the stopping time by \E[\E[K | \rho=z]^{\gamma}], where K is a random stopping time
+    stop_time_log_part = np.sqrt(1 + gamma)
+    stop_time_log_part = stop_time_log_part * np.power(np.sqrt(2 * np.pi), gamma)
+    stop_time_log_part = stop_time_log_part * ( np.power(( (margin * (1 + gamma))/sigma_1), gamma ) + 1)
+    stop_time_log_part = np.log(stop_time_log_part)
+    stop_time_exp_part = (gamma * margin**2) / (2 * sigma_1**2)
+
+    # note exp(0) = 1, which is implicitly added
+    moment_term = utils.stable_logsumexp_two(0, stop_time_log_part + stop_time_exp_part)
+    moment_term = moment_term / (gamma * (alpha - 1))
+
+    moment_based_bound = rdp_rho + rdp_nu + moment_term
+    return min(moment_based_bound, ret_rdp)
 
 
 

From 6df642f0c31ce967ee74ac4fe412eeba3e0d165d Mon Sep 17 00:00:00 2001
From: Jonathan Lebensold <lebensold@meta.com>
Date: Fri, 7 Jul 2023 15:22:36 -0400
Subject: [PATCH 3/8] revert to master

---
 autodp/rdp_bank.py                     | 38 ++++++++++++++++----------
 papers/NeurIPS20-SVT/exp2_align_var.py | 19 ++++++-------
 2 files changed, 31 insertions(+), 26 deletions(-)

diff --git a/autodp/rdp_bank.py b/autodp/rdp_bank.py
index f81e107..7e9f241 100644
--- a/autodp/rdp_bank.py
+++ b/autodp/rdp_bank.py
@@ -6,13 +6,19 @@
 
 Some of the functions contain the renyi divergence of two given distributions, these are useful to keep track of
 the per-instance RDP associated with two given data sets.
+
+TO CONTRIBUTORS:  1. any new addition to the rdp_bank should include a reference to the mechanism of
+    interesting and the derivation of its RDP.
+    2. You should try providing an implementation fo the entire range of alpha >0. If you do not
+       have alpha <1, feel free to use any upper bound, e.g., the bound for alpha = 1.
+
 """
 
 
 import numpy as np
 import math
 from autodp import utils
-from autodp.utils import _log1mexp, logcomb, stable_logsumexp
+from autodp.utils import _log1mexp
 
 
 def stable_log_diff_exp(x):
@@ -50,13 +56,20 @@ def RDP_laplace(params, alpha):
     # assert(b > 0)
     # assert(alpha >= 0)
     alpha=1.0*alpha
-    if alpha <= 1:
-        return (1 / b + np.exp(-1 / b) - 1)
-    elif np.isinf(alpha):
+
+    if np.isinf(alpha):
         return 1/b
-    else:  # alpha > 1
+    elif alpha == 1:
+        # KL-divergence
+        return 1 / b + np.exp(-1 / b) - 1
+    elif alpha > 1:  # alpha > 1
         return utils.stable_logsumexp_two((alpha-1.0) / b + np.log(alpha / (2.0 * alpha - 1)),
                                            -1.0*alpha / b + np.log((alpha-1.0) / (2.0 * alpha - 1)))/(alpha-1)
+    elif alpha == 0.5:
+        return -2*(-1.0/(2*b) + np.log(1 + 1.0/(2*b)))#   -2*np.log(np.exp(-1.0/(2*b))*(1+1.0/(2*b)))
+    else:
+        return np.log(alpha/(2.0*alpha-1)*np.exp((alpha-1.0)/b) + (alpha-1.0)/(2.0*alpha-1)*np.exp(-1.0*alpha/b))/(alpha-1)
+        # Handling the case when alpha = 1/2?
 
 def RDP_zCDP(params,alpha):
     """
@@ -359,20 +372,15 @@ def RDP_gaussian_svt_cgreater1(params, alpha):
 
     Args:
         k:the maximum length before svt stops
-        sigma: noise added to the threshold divide by L2 sensitivity. 
-		sigma_nu: noise added to the query divide by 2 * L2 sensitivity.
+        sigma: noise added to the threshold.
         c: the cut-off parameter in SVT.
     """
-    sigma_rho = params['sigma']
-    sigma_nu = params['sigma_nu']
+    sigma = params['sigma']
     c = max(params['c'], 1)
     k = params['k']  # the algorithm stops either k is achieved or c is achieved
-    rdp_rho = 0.5 / (sigma_rho ** 2) * alpha
-    rdp_nu = 0.5 / (sigma_nu ** 2) * alpha
-    #c_log_n_c = c * np.log(k / c) # approximate bound
-    log_comb = [logcomb(k, i) for i in range(c+1)]
-    log_sum_comb = stable_logsumexp(log_comb)
-    ret_rdp = log_sum_comb * 1.0 / (alpha - 1) + rdp_rho + rdp_nu * c
+    rdp_rho = 0.5 / (sigma ** 2) * alpha
+    c_log_n_c = c * np.log(k / c)
+    ret_rdp = c_log_n_c * 1.0 / (alpha - 1) + rdp_rho * (c + 1)
     return ret_rdp
 
 
diff --git a/papers/NeurIPS20-SVT/exp2_align_var.py b/papers/NeurIPS20-SVT/exp2_align_var.py
index ddd6cb5..cd17a3b 100644
--- a/papers/NeurIPS20-SVT/exp2_align_var.py
+++ b/papers/NeurIPS20-SVT/exp2_align_var.py
@@ -7,7 +7,7 @@
 import matplotlib.pyplot as plt
 from scipy.stats import norm, laplace
 from scipy.special import comb
-#import matplotlib.font_manager as fm
+import matplotlib.font_manager as fm
 from autodp.mechanism_zoo import ExactGaussianMechanism, PureDP_Mechanism,SubsampleGaussianMechanism, GaussianMechanism, ComposedGaussianMechanism,GaussianSVT_Mechanism, NoisyScreenMechanism
 from autodp.transformer_zoo import Composition, AmplificationBySampling
 
@@ -26,7 +26,6 @@
 delta = 1e-6
 lambda_rho = 120
 lambda_nu = 240
-# Under this setting, sigma_1 / L2 sensitivity = sigma_2 /2 * L2 sensitivity
 sigma_1 = lambda_rho*np.sqrt(2)
 sigma_2 = 2*sigma_1
 eps_1 = 1.0 / lambda_rho
@@ -103,12 +102,11 @@ def exp_2a():
 
             # stage-wise generalized SVT, k is the maximum length of each chunk
             k = int(idx / np.sqrt(count_gau))
-            # assume sensitivity is 1, sigma_nu denotes noise add to query / 2
-            generalized_mech = StageWiseMechanism({'sigma':sigma_1,'k':k, 'c':count_gau, 'delta':delta})
+            generalized_mech = StageWiseMechanism({'sigma':sigma_1,'k':k, 'c':count_gau})
             eps_kov.append(generalized_mech.get_approxDP(delta))
 
             # Gaussian-SVT c>1 with RDP, k is the total length before algorithm stops
-            gaussianSVT_c = GaussianSVT_Mechanism({'sigma':sigma_1,'sigma_nu':sigma_1,'k':idx, 'c':count_gau}, rdp_c_1=False)
+            gaussianSVT_c = GaussianSVT_Mechanism({'sigma':sigma_1,'k':idx, 'c':count_gau}, rdp_c_1=False)
             eps_g.append(gaussianSVT_c.get_approxDP(delta))
 
             #Gaussian-SVT with c=1, we use average_k as the approximate maximum length of each chunk, margin is used in Proposition 10
@@ -116,7 +114,6 @@ def exp_2a():
             params_SVT = {}
             params_SVT['k'] = average_k
             params_SVT['sigma'] = sigma_1
-            params_SVT['sigma_nu'] = sigma_1
             params_SVT['margin'] = margin
             per_gaussianSVT_mech = GaussianSVT_Mechanism(params_SVT)
             gaussianSVT_mech = compose([per_gaussianSVT_mech],[max(count_gau, 1)])
@@ -132,7 +129,7 @@ def exp_2a():
             'weight': 'bold',
             'size': 18}
 
-    #props = fm.FontProperties(family='Gill Sans', fname='/Library/Fonts/GillSans.ttc')
+    props = fm.FontProperties(family='Gill Sans', fname='/Library/Fonts/GillSans.ttc')
     f, ax = plt.subplots()
     plt.figure(num=0, figsize=(12, 8), dpi=80, facecolor='w', edgecolor='k')
     plt.loglog(index, eps_a, '-r', linewidth=2)
@@ -151,7 +148,7 @@ def exp_2a():
     plt.yticks(fontsize=20)
     plt.xlabel(r'Iterations', fontsize=20)
     plt.ylabel(r'$\epsilon$', fontsize=20)
-    #ax.set_title('Title', fontproperties=props)
-    plt.show()
-    #plt.savefig('exp2a.pdf', bbox_inches='tight')
-exp_2a()
+    ax.set_title('Title', fontproperties=props)
+
+    plt.savefig('exp2a.pdf', bbox_inches='tight')
+exp_2a()
\ No newline at end of file

From 68f09a23e3491800ed49904d7c3a76b0380ced8b Mon Sep 17 00:00:00 2001
From: Jonathan Lebensold <lebensold@meta.com>
Date: Fri, 7 Jul 2023 15:26:02 -0400
Subject: [PATCH 4/8] match master

---
 autodp/dp_bank.py       | 19 ++++---------------
 autodp/mechanism_zoo.py |  4 ++--
 2 files changed, 6 insertions(+), 17 deletions(-)

diff --git a/autodp/dp_bank.py b/autodp/dp_bank.py
index bcc0028..1af1456 100644
--- a/autodp/dp_bank.py
+++ b/autodp/dp_bank.py
@@ -85,33 +85,22 @@ def fun(x):
     else:
         raise RuntimeError(f"Failed to find epsilon: {results.flag}")
 
-def eps_generalized_gaussian(params, delta):
+def eps_generalized_gaussian(x, sigma, delta,k, c, c_tilde):
     """
     submodule for generalized SVT with Gaussian noise
     we want to partition c into [c/c'] parts, each part using (k choose c')
     need to check whether (k choose c') > log(1/delta')
     k is the maximam number of queries to answer for each chunk
+    x is log delta for each chunk, it needs to be negative
     :param x:
     :param sigma:
     :param delta:
     :return:
     """
-    sigma = params['sigma']
-    k = params['k']
-    c = params['c']
-    # delta' (delta1) is the delta budget for each chunk.
-    c_tilde = int(np.sqrt(c))
-    # divide half of the delta budget into c_tilde chunks, x is log delta for each chunk, it needs to be negative.
-    x = np.log(delta/(2 * c_tilde))
     acct = dp_acct.DP_acct()
+    per_delta = np.exp(x) # per_delta for each c' chunk
     coeff = comb(k,c_tilde)
-    # Find a proper per delta such that c_tilde * per_delta < delta and (k choose c') > log(1/per_delta)
-    per_delta = min( 0.5 /coeff, delta/(2*c_tilde))
-    x = np.log(per_delta)
-    if per_delta > 1.0/(coeff):
-        print('per delta', per_delta, 'coeff', coeff, 'k', k, 'c', c)
     assert per_delta < 1.0/(coeff)
-
     #compute the eps per step with 1/(sigma_1**2) + sqrt(2/simga_1**2 *(log k + log(1/epr_delta)))
     # compose eps for each chunk
     while c:
@@ -156,4 +145,4 @@ def get_eps_randresp_optimal(p,delta):
     elif delta >= 2*p - 1:
         return 0.0
     else:
-        return np.log(p-delta) - np.log(1 - p)
+        return np.log(p-delta) - np.log(1 - p)
\ No newline at end of file
diff --git a/autodp/mechanism_zoo.py b/autodp/mechanism_zoo.py
index 12a00ec..ce7ca1d 100644
--- a/autodp/mechanism_zoo.py
+++ b/autodp/mechanism_zoo.py
@@ -385,11 +385,11 @@ def __init__(self, params=None,approxDP_off=False, name='StageWiseMechanism'):
         Mechanism.__init__(self)
 
         self.name = name # When composing
-        self.params = params
+        self.params = {'sigma': params['sigma'], 'k':params['k'], 'c':params['c']}
         self.delta0 = 0
 
         if not approxDP_off:  # Direct implementation of approxDP
-            new_approxdp = lambda x: dp_bank.eps_generalized_gaussian(params, x)
+            new_approxdp = lambda x: dp_bank.get_generalized_gaussian(params, x)
             self.propagate_updates(new_approxdp, 'approxDP_func')
 
 

From b8f1fcd1618e2df38cc574dfd3f023d37e53e29c Mon Sep 17 00:00:00 2001
From: Jonathan Lebensold <lebensold@meta.com>
Date: Fri, 7 Jul 2023 15:29:44 -0400
Subject: [PATCH 5/8] match master

---
 autodp/mechanism_zoo.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/autodp/mechanism_zoo.py b/autodp/mechanism_zoo.py
index ce7ca1d..cb24bc9 100644
--- a/autodp/mechanism_zoo.py
+++ b/autodp/mechanism_zoo.py
@@ -352,7 +352,7 @@ def __init__(self,params,name='GaussianSVT', rdp_c_1=True):
             new_rdp = lambda x: rdp_bank.RDP_gaussian_svt_c1(self.params, x)
         else:
             self.name = name + 'c>1'
-            self.params = {'sigma':params['sigma'],'sigma_nu': params['sigma_nu'], 'k':params['k'], 'c':params['c']}
+            self.params = {'sigma':params['sigma'],'k':params['k'], 'c':params['c']}
             new_rdp = lambda x: rdp_bank.RDP_gaussian_svt_cgreater1(self.params, x)
         self.propagate_updates(new_rdp, 'RDP')
 

From cf1e200ab0b1e875394e460e67f8d9768080b477 Mon Sep 17 00:00:00 2001
From: Jonathan Lebensold <lebensold@meta.com>
Date: Fri, 7 Jul 2023 15:31:21 -0400
Subject: [PATCH 6/8] whitespace

---
 autodp/dp_bank.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/autodp/dp_bank.py b/autodp/dp_bank.py
index 1af1456..b04afaa 100644
--- a/autodp/dp_bank.py
+++ b/autodp/dp_bank.py
@@ -145,4 +145,4 @@ def get_eps_randresp_optimal(p,delta):
     elif delta >= 2*p - 1:
         return 0.0
     else:
-        return np.log(p-delta) - np.log(1 - p)
\ No newline at end of file
+        return np.log(p-delta) - np.log(1 - p) 
\ No newline at end of file

From 2d787a0d722204a95c969ab015dba42e90e62e97 Mon Sep 17 00:00:00 2001
From: Jonathan Lebensold <lebensold@meta.com>
Date: Fri, 7 Jul 2023 15:32:11 -0400
Subject: [PATCH 7/8] whitespace

---
 autodp/dp_bank.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/autodp/dp_bank.py b/autodp/dp_bank.py
index b04afaa..1af1456 100644
--- a/autodp/dp_bank.py
+++ b/autodp/dp_bank.py
@@ -145,4 +145,4 @@ def get_eps_randresp_optimal(p,delta):
     elif delta >= 2*p - 1:
         return 0.0
     else:
-        return np.log(p-delta) - np.log(1 - p) 
\ No newline at end of file
+        return np.log(p-delta) - np.log(1 - p)
\ No newline at end of file

From 972dc84b977ce36c87c6c3b6be68a19d4fc1c97b Mon Sep 17 00:00:00 2001
From: Jonathan Lebensold <lebensold@meta.com>
Date: Fri, 7 Jul 2023 15:33:54 -0400
Subject: [PATCH 8/8] whitespace

---
 autodp/dp_bank.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/autodp/dp_bank.py b/autodp/dp_bank.py
index 1af1456..4fbcbaf 100644
--- a/autodp/dp_bank.py
+++ b/autodp/dp_bank.py
@@ -145,4 +145,4 @@ def get_eps_randresp_optimal(p,delta):
     elif delta >= 2*p - 1:
         return 0.0
     else:
-        return np.log(p-delta) - np.log(1 - p)
\ No newline at end of file
+        return np.log(p-delta) - np.log(1 - p)