Diffusion Sampling with Momentum for Mitigating Divergence Artifacts
Authors: Suttisak Wizadwongsa
Worameth Chinchuthakun
Pramook Khungurn
Amit Raj
Supasorn Suwajanakorn
Venue: ICLR 2024 poster (2024)
Keywords: diffusion models
heavy ball momentum
divergence artifacts
numerical method
ode solver
image generation
Primary area: generative models
TL;DR: This paper addresses the issue of unexpected divergence artifacts in diffusion sampling and proposes two novel techniques to address them

Abstract: Despite the remarkable success of diffusion models in image generation, slow sampling remains a persistent issue. To accelerate the sampling process, prior studies have reformulated diffusion sampling as an ODE/SDE and introduced higher-order numerical methods. However, these methods often produce divergence artifacts, especially with a low number of sampling steps, which limits the achievable acceleration. In this paper, we investigate the potential causes of these artifacts and suggest that the small stability regions of these methods could be the principal cause. To address this issue, we propose two novel techniques. The first technique involves the incorporation of Heavy Ball (HB) momentum, a well-known technique for improving optimization, into existing diffusion numerical methods to expand their stability regions.  We also prove that the resulting methods have first-order convergence. The second technique, called Generalized Heavy Ball (GHVB), constructs a new high-order method that offers a variable trade-off between accuracy and artifact suppression. Experimental results show that our techniques are highly effective in reducing artifacts and improving image quality, surpassing state-of-the-art diffusion solvers on both pixel-based and latent-based diffusion models for low-step sampling. Our research provides novel insights into the design of numerical methods for future diffusion work.

Published as a conference paper at ICLR 2024
DIFFUSION SAMPLING WITH MOMENTUM FOR MITI- GATING DIVERGENCE ARTIFACTS
Suttisak Wizadwongsa
Worameth Chinchuthakun Tokyo Institute of Techhnlogy,, Japan ch mrinnhhthaaeee hwnuuthaaunn... w.a.m. titech. ac. ip
VISTEC, Thailand
suttisak. s1 1.9.vitte........
Pramook Khungurn pixiv Inc. pr zramook@gmaile com
Amit Raj Google
Supasorn Suwajanakorn
VISTEC, Thailand supasorn. s@vistec.
amitrajs@@ooglee com
ABSTRACT
Despite the remarkable success of diffusion models image generation, slow sampling remains a persistent issue. To accelerate the sampling process, prior studies have reformulated diffusion sampling as an ODE/SDE and introduced higher-order numerical methods. However, these methods often produce diver- gence artifacts, especially with a low number of sampling steps, which limits the achievable acceleration. In this paper, we investigate the potential causes of these artifacts and suggest that the small stability regions of these methods could be the principal cause. To address this issue, propose two novel techniques. The first technique involves the incorporation Heavy Ball (HB) momentum, a well-known technique for improving optimization, into existing diffusion numer- ical methods to expand their stability regions. We also prove that the resulting methods have first-order convergence. The second technique, called Generalized Heavy Ball (GHVB), constructs a new high-order method that offers variable trade-off between accuracy and artifact suppression. Experimental results show that our techniques are highly effective in reducing artifacts and improving im- age quality, surpassing state-of-the-art diffusion solvers on both pixel- pixel-based and latent-based diffusion models for low-step sampling. Our research provides novel insights into the design of numerical methods for future diffusion work.
INTRODUCTION
Diffusion models (Ho al. 2020 Song al. 2020a) family of generative models that has gar- nered considerable attention due to their remarkable image quality. Unlike Generative Adversarial Networks (GANs) (Goodfellow et al., 2014), which may suffer from mode collapse and instabilities
during training, diffusion models are less sensitive to hyperparameters (Ho et al.] 22200 Rombach al. 2022) and offer better sampling quality (Dhariwal & Nichol 2021). These models have been successfully applied to various image-related tasks, such as text-to-image generation (Nichol et al.)
2022). image-to-image translation (Su al. (002), image composition (Sasaki al. [021]. adver-
sarial purification (Wang et al., 2022] Wu et al. 2022), and super-resolune resouutisoo (Choi et al. (2021)
However, one significant drawback of diffusion models their slow sampling speed. This is be- cause the sampling process involves Markov chain, which requires a large number of iterations to generate high-quality results. Recent attempts to accelerate this process include improvements the noise schedule (Nichol & Dhariwal 2021 Watson et al.) 202T) and network distillation (Salimans) & Ho 2022 Watson et 2022 Song et al. 2023). Fortunately, the sampling process can be repre- sented by ordinary stochastic difterential equations, and numerical methods can be used to reduce the number of iterations required. While DDIM (Song et al., 2020a)), a 1 method, the most commonly used approach, it still requires considerable number of iterations. Higher-orde numer- ical methods, such as DEIS (Zhang & Chen 2022). DPM-Solver (Lu et al.] (202a)), and PLMS (Liu al. 2022a), have been proposed generate high-quality images in fewer steps. However, these

---

Published as a conference paper at ICLR 2024
DPM-Solver++ (2M) Lu al.
2022a
PLMS4 et al.
2022a
Without Momentum
With HBB 0.5 (Ours)
GHVB 1.5 (Ours)
15
25
20
Figure 1: We demonstrate the occurrence of divergence artifacts in existing numerical methods at various guidance scales & using 15 sampling steps. Integrating HB momentum into these methods and our GHVB 1.5 method can mitigate these artifacts. Prompt: "photo of girl face" (a) Realistic Vision v2.0, (b) Anything Diffusion v4.0, (c) Deliberate Diffusion.
methods begin to produce artifacts (see Figure D when the number of steps is decreased beyond a certain value, thereby limiting how much we can reduce the sampling time.
In this study, we investigate the potential causes of these artifacts and found that the narrow stabil- ity region of high-order numerical methods cause solutions to diverge, resulting in divergence artifacts. To address this issue and enable low-step, artifact-free sampling, we propose two tech- niques. The first technique involves incorporating Heavy Ball (HB) momentum (Polyak) 1987), well-known technique for improving optimization, into existing diffusion numerical methods. This approach effectively reduces divergence artifacts, but its accuracy only has first order of conver- gence. In this context, the accuracy measures how close the approximated, low-step solution is to the solution computed from very high-step solver (e.g.. 1,000-step DDIM). The second technique, called Generalized Heavy Ball (GHVB), is new high-order numerical method that offers vari- able trade-off between accuracy and artifact suppression. Both techniques are training-free and ineur negligible additional computational costs. [igure[[demonstrrr reducing artifacts compared to previous diffusion sampling methods. Furthermore, our experiments show that our techniques are effective both pixel-basee and latent-based diffusion models.
the capability of both techniques in
BACKGROUND
This section first presents the theoretical foundation of diffusion sampling when modeled an ordinary differential equation (ODE) and related numerical methods. Second, we discuss ODE forms for guided diffusion sampling and prior splitting numerical methods. Third, we review the concept of stability region, which our primary analysis tool.
2.1 DIFFUSION IN ODE FORM
Modeling diffusion sampling as an ODE commonly based on the non-Markovian sampling of Denoising Diffusion Implicit Model (DDIM) from Song et al. (2020a). DDIM is given by:
04--1 (xr /1 (x,o( + V a1-160(111
(1)
Here, (0 t) is neural network with learnable parameters 0 trained to predict noise given the current state x, and time as input. The parameter a, is a schedule that controls the degree of diffusion at each time step. Previous research has shown that DDIM (T) can be rewritten into
2

---

Published as a conference paper at ICLR 2024
ODE, making it possible to use numerical methods to accelerate the sampling process. Two ODEs have been proposed in the literature:
di
dē O( dē
=
(2)
do
Equation 2 can be obtained by pe-prammteriiingg 0 =V- at. I zz pr/van and a(F 0)
<g( (). These transformations are widely used in various diffusion solvers (Zhang & Chen 2022 Lu 2022a Liu al. 2022a Zhao et al. 2023). If (are, 1) is sum of multiple terms, such as in split Equation and solve each resulting equation separately (Wizadwongsal & Suwajanakorn 2023). Another ODE (Equation 3) is derived by defining = var/ // a and
guided diffusion,
21/ 21/v1 at where an estimation of the final result s(ē, ō) = (æ; v1 - axc t))/ vai- This ODE has the advantage of keeping the differentiation bounded within the pixel value range in pixel-based diffusion, outperforning Equation 2 in many cases in DPM-Solver++ (Lu et al. 2022b
2.2 STABILITY REGION
The stability region is fundamental concept in numerical methods for solving ODEs. determines the step sizes that enable numerical approximations to converge. To illustrate this concept, let us consider the Euler method, simple, ordorder method for solving ODEs, given by
Tn+11 = In + of(rn),
(4)
where n, is the approximate solution and 5 is the step size. To analyze the stability of the Euler method, we can consider a test equation of the form x' = Ax, where is a complex constant. The solution of this test equation can be expressed
In+1 = In + бλτ, = (1 + 8A)xn, = (1 + 8A)n+1ao 10.
(5)
where Zo is the initial value. For the approximate solution to converge to the true solution, it is necessary that | + SA I 1. Hence, the stability region of the Fuler method S (z E C [1 + =| & }} because if 5A lies outside of S, the solution In will tend to too as n DO.
In diffusion sampling, another common family of numerical solvers is the to as Pseudo Linear Multi-Step (PLMS). AB methods en- compass the Fuler method as its 1storder (AB1), and the 2nd -oderr (AB2) is given by:
(AB) methods, also referred
1.0
AB1 (Eulor)
AB3
A84
3 2
1 3 - )
En+11 Ir, 6
-05-
The stability regions of AB methods of various orders are derived in Appendix C These regions can be visualized by determining their boundaries using the boundary locus tech-
-1.0
33.
-25
-2.0
-1.5
-1.0
-0.5
0.0
nique (Lambert et al. [991). As depicted in Figure [ the stability region decreases in size and its boundary becomes more restrictive as the order of the method increases.
Figure 2: Boundaries of stability re- gions of the first Adams-Bashtone methods.
UNDERSTANDING ARTIFACTS IN DIFFUSION SAMPLING
One unique issue in diffusion sampling is the occurrence of "divergence" artifacts, which are char- acterized by regions with unrealistic, oversaturated pixels in the outputs. This problem can arise from several factors including the use of high-order numerical solvers, too few sampling steps, or high guidance scale. (See Appendix (P). The existing solution is simply avoiding these factors at the cost of slower sampling speed or weaker guidance. This section investigates the cause of these artifacts in Section 3.1 and proposes solutions that do not sacrifice sampling speed in Section 3.2
3.1 ANALYZING DIFFUSION ARTIFACTS
We analyze the areas where divergence artifacts occur during sampling by examining the magnitudes of the latent variables in those areas. Specifically, we Stable Diffusion (Rombach et al. 2022
3

---

Published as a conference paper at ICLR 2024
PLMS4 250 steps
(b) PLMS4 steps
Figure 3: Generated images and max- poxxepooed latent magnitude maps (brighter colors indicate higher values). Images with artifacts (b) have higher latent magnitudes than those without artifacts (a).
which operates and performs diffusion sampling on a latent space of dimension 64 x 64 x 4, to generate images with and without artifacts by varying the number of steps. Then, visualize each latent variable E RRd in the 64 x 64 spatial grid after normalizing it using the channel-wise mean and standard deviation, computed from the COCO dataset (Lin et al.) 2014). Figure shows the magnitudes of the normalized latent variables after 4 X 4 max pooling for visualization purposes.
We found that artifacts mainly appear in areas where the latent magnitudes are higher than usual. Note that images without artifacts can also have high latent magnitudes in some regions, although this is exceptionally rare. Conversely, when artifacts appear, those regions almost always have high magnitudes. In pixel-based diffusion models, the artifacts manifest directly as pixel values near I or due to clipping, which can be observed in Figure 17 in Appendix
3.2 CONNECTION BETWEEN ODE SOLVER AND ARTIFACTS
We hypothesize that these artifacts are caused by the numerical instability during sampling. Math- ematically, analyze the ODE for diffusion sampling in Equation 2 using the problem reduction technique for stiffness analysis (Higham & Trefethen 1993). Assuming that the effect of a on the function is negligible, use Taylor expansion approximate the RHS. of Equation yielding
dæ do
Ve(=") x²) O(IE x 1²1)²).
Here, x" denotes the converged solution that should not have any noise left (i.c. E(x") = 0), and Vε(x*) denotes the Jacobian matrix at x². As converges x², the term O[ x* 2=| becomes negligibly small, so we may drop it from the equation.
Let À be an eigenvalue of Vē(=*)² and v be the corresponding normalized eigenvector such that Vē(x*)" Av. We define u = = (x - -x") and obtain n' Au as our test equation. According Section if 8A falls outside the stability region of numerical method, the numerical solution to may diverge, resulting in diffusion sampling results with larger magnitudes that later manifest as divergence artifacts. Therefore, divergence artifacts are more likely to occur when the stability region is too small. Although some numerical methods have infinite stability regions, those used in diffusion sampling have only finite stability regions, which implies that the solution will always diverge if the step size 8 is sufficiently high. More details about the derivation Appendix D and a 2D toy example illustrating this effect is provided in Appendix E
be found in
One possible solution to mitigate artifacts is to reduce the step size 8, which shifts 8A closer to the origin of the complex plane. However, this approach increases the number of steps, making the process slower. Instead, we will modify the numerical methods to enlarge their stability regions.
4

---

Published as a conference paper at ICLR 2024
4 METHODOLOGY
This section describes two techniques for improving stability region and reducing divergence arti- facts. Specifically, first show how to apply Polyak's Heavy Ball Momentum (HB) to diffusion sampling, and secondly, how to generalize it to higher orders. Our techniques are designed to be simple to implement and do not require additional training.
4.1 POLYAK'S S HEAVY BALL MOMENTUM FOR DIFFUSION SAMPLING
Recall that Heavy Ball Momentum is an optimization algorithm proposed by Polyak (1987) that enhances gradient descent (In+1 2n 3, f )). The method takes inspiration from the physical analogy of a heavy ball moving through a field of potential with damping friction. The update rule for Polyak's HB optimization algorithm is given by:
Fn+1 = In + Onn ( In-1) 3, (x,))
(8)
where On and 3n are parameters. We can apply HB to the Euler method in which case we typically set an (1 (1,) obtain
Tn+1 = In + 3,) - xn-1) - (z.))
and we may show that the numerical method above has the same order of convergence as the original Euler method. For simplicity, we assume that 5, 5 E (0, 1), which is constant known as the damping coefficient. Then, we can reformulate Equation as:
Un+1 (1 3) + (f()),
En+1 In dUn+1
(10) the exponential moving average
Here, we may interpret v, as the heavy ball's position, and of ( its velocity. We can see that position is updated with "displacement time velocity," much like in physics.
For high-order method of the form Int1 = =In + Ek (x we can apply HB to it as:
Un+1 (1 B)e + 3 Σ bif
Pnt1 In + ōun+1.
(11)
1=0
The resulting numerical method has larger stability region, as can be seen in Figures 4b to 4d in
which we show stability boundaries of methods after HB is applied to them with varying fs.. (We use HB 0.4 to denote 3 0.4). However, Theorem in Appendix shows that soon as deviates from 1, the theoretical order of convergence drops to I, leading to a significant decrease in image quality, illustrated in F In the next subsection, we propose an alternative approach that increases the stability region while maintaining high order of convergence.
65-
-20
-2.
-10
-20
-25
-20
-35
-1a
-55
-20
-25
-20
23
-22
(a) order
(b) 2dd order
(c) 3 order
(d) order
(PL.M w/ HR
04-
(PI MS2 HR 0.41.)
(PI M53 w/ 0.- L.)
(PI MS4 w/ HR 0.4
Figure 4: Boundaries stability regions of 1s order AB methods with HB applied to them with different values of the damping coefficient B.
4.2 GENERALIZING POLYAK'S ' HEAVY BALL TO HIGHER ORDERS
In this section, generalize the Euler method with HB momentum achieve high-order conver- gence, similar to AB methods' generalization of the Fuler method. Defining the backward difference
operator A Az, I, In-1- Following Berry & Healy (2004). we express the AB formula by:
3
251
95
Ar 5 1 A 2
A
A +...) ((n))
(12)
720
288
5

---

Published as a conference paper at ICLR 2024
B 0.2
= 0.1
B 0.6
0.8
PLMSI (3 1.0)
(a) PLM54 with HB
B 0.2
= 0.4
, 0.6
0.8
I3D (1,000 stops)
(b) 4thh order GHVB
Figure 5: Comparison between the two proposed techniques: (a) HB and (b) GHVB, applied to PLMS4 (Liu al., 2022a) with 15 sampling steps. Both effective at reducing artifacts, but HB's ' accuracy drops faster than GHVB's ' as 3 deviates from 1. Positions of the lanterns in Row (a) deviate more from the ground-truth (1,000-steps DDIM (Zhang et al. 2023)) than those in Row (b). Moreover, the image at 0.2 in Row (a) becomes blurry as HB yields numerical method with a lower order of convergence than what GHVB does. Prompt: "A beautiful illustration of people releasing lanterns near a river".
The order convergence is determined by the number of terms on the RHS; for example, AB2 can be written as Ar (1 + 1 f(ra). Multiplying Equation 12 by (B + (1 3)4), we have:
2 3
6
10
(3+ + 3) =
+ .... (xx)).
(13)
2
12
24
Note that the update of on in Equation 10 can be rewritten as (3+(1 (1- 3) = =Bf(z)) Next, choose the order of convergence by fixing the number of terms the RHS. To get, say, -od-rr method, we may choose:
2
B
2 A 3f(xn)
(B + (1 -1-11 ()-3)))Aam.. =
B
A
[(xn) =
(14)
2
23
2
=8 1
A
(3 + (1 (1-8)A)en+1.
(15)
23
Eliminating (B - -)) from both sides, we obtain the 2nd-order generalized HB method:
Mn+1 (1 8)v, + Bf(xn), Int1 + 5 ( 2 + Unt1 +
2
(16)
2,9
23
Algorithm details a complete implementation. When 3 = 1, the formulation in Equation 16 is equivalent to the AB2 formulation in Equation 6 As approaches O, Equation [4 converges the 1--order Euler method 4 Thus, this generalization also serves as an interpolating technique between two adjacent-order AB methods, except for the ordorder GHVB, which equivalent to the Euler method with HB momentum in Equation 10
We call this new method the Generalized Heavy Ball (GHVB) and associate with it momentum number, whose ceiling indicates the method's order. For example, GHVB 1.8 refers to the 2nd -order GHVB with 3 = 0.8. The main difference between HB and GHVB is that HB caleulates the moving average after summing high-order coefficients, whereas GHVB calculates it before the summation.
We analyze the stability region of GHVB using the same approach as before and visualize the re- gion's 's locus curve in Figure The theoretical order of accuracy of this method is given by Theorem in Appendix(H) [ We discuss alternative momentum methods, such as Nesterov's momentum, which can offer comparable performance but are less simple in Appendix G
5 EXPERIMENTS
We present a series of experiments to evaluate the effectiveness of our techniques. In Section 5.1 we assess the divergence artifacts reduction through qualitative results and quantitative measurements

---

Published as a conference paper at ICLR 2024
65
63
65
03-
80
25
20
14
5.3
-20
-25
-10
-1
-10
-35
10
-0
-2
-20
455
-0
-05
60
-30
-15
-21
-35
80
(a) 1d order
(b) 2 order
(c) 3dd order
(d) 4 order
(GHVB 0.-
(GHVB 1.4
(GHVB 244
(GHVB 3.4 4.00
Figure 6: Boundary of stability regions for 1 4th order Generalized Heavy Ball (GHVB).
of the latent magnitudes in text-to-image diffusion model. Besides reducing artifacts, another important goal is to ensure that the overall sampling quality improves and does degenerate (e.g., becoming color blobs). We test this with experiments on both pixel-based and bheent diffusion models trained on ImageNet256 (Russakovsky et al. 2015) (Section 5.2 and (.3). which show that our techniques significantly improve image quality, as measured by the standard Fréchet Inception Distance (FID) score. Lastly, Section 5.4] we present an ablation study of GHVB methods with
varying degrees of order. A similar study on HB methods can be found in Appendix Q
3.0
PLM54 HB D.3
PLAIS
3.0
PLMS5
GHVB37
DAN Solvere+ wi HB
Solvare+
2.5
GHVB2.7
GIVR333
USSP (PLMS4
DPM-Soleere+ DPM Solvort- Me1 0.8
FD Score
2.0
2.0
GHVVBL.
DOIN
60
3.5
1.2
10
20
40
50
10
30
30
60
10
15
22
26
Number Staps
Number Steps
Number of Stops
Figure 7: Mean magnitude scores
Figure 8: Mean magnitude scores Figure 9: FID on ADM (lours)
5.1 ARTIFACTS MITIGATION
In this experiment, we apply our HB and GHVB techniques to the most popular 2nd and 4th -oddrr solvers, DPM-Solver++ (Lu et al., (022a) and PLMS4 (Liu et al. (2022a). using 15 sampling steps and various guidance scales on three different text-to-image diffusion models. The qualitative re- sults in Figure show that our techniques significantly reduce the divergence artifacts and produce realistic results (columns a, c). More qualitative results are in Figure 16 in Appendix! J
Quantitatively measuring divergence artifacts can be challenging, metrics like MSE LPIPS may only capture the discrepancy between the approximated and the true solutions, which does not necessarily indicate the presence of divergence artifacts. In this study, use the magnitudes of latent variables as introduced in Section 3.1 as proxy metric to measure artifacts. In particular, we define magnitude score
Ei, ( that sums over the latent variables in max max-pooled latent grid, where g(x) x if 2 7 and 0 otherwise. We generate 160 samples from the same set of text prompts and seeds for each method from fine-tuned Stable Diffusion model called Anything V4.
The results using 7 3 (magnitude considered high when above std.) are shown in Figures 7 and We observe that the magnitude score increases as the number of sampling steps decreases and orgeer-order methods result in higher magnitude scores. Adding HB momentum PLMS4 (Liu et al. 2022a) DPM-Solver++ (Lu et al., (022b) (Figure or employing GHVB with momentum number 3 < (Figure can reduce these latent magnitudes. Our methods also show larger per- formance difference from previous works the problem setup is more challenging, as demonstrated in additional Stable Diffusion experiments in Appendix M and
8
Our recommendation is to employ the Euler method with HB 0.8 or GHVB 0.8 the default choices for common usage due to their consistent performance and sufficiently expansive stability regions (see Appendix for more experiments). In scenarios where accuracy fidelity to the true DDIM solution is essential or when lower sampling step is required, higher-order GHVB variants can be considered. For guidance method selection, please refer Q2 in Appendix

---

Published as a conference paper at ICLR 2024
5.2 EXPERIMENTS ON PIXEL-BASED DIFFUSION MODELS
We evaluate techniques using Klasifier-guiiee diffusion sampling with ADM from Peebles & Xie (2022), an unconditioned pixel-based diffusion model, with their classifier model. Additionally, compare methods with other diffusion sampling methods, namely DPM-Solver++ (Lu et al. (022b) and LTSP (Wizadwongsa & Suwajanakorn, 2023), which have demonstrated strong performance in classifier-guided sampling.
For DPM-Solver++, we use -odder multi-step method and compare the results with and without HB momentum. For LTSP, split numerical method, use PLMS4 (Liu al.) 2022a) to solve the first subproblem (see (Wizadwongsa & Suwajanakorn, 2023)) and compare different methods for solving the second subproblem, including regular Euler method and its variant with HB momentum (equivalent to GHVB 0.8). Figure Indiiates that our techniques effectively improve FID scores for both + and LTSP. Notably, applying our HB momentum to LTSP consistently yields the lowest FID scores. This experiment highlights the benefits of using HB momentum, which presents additional results, and Figure
provides better choice than Euler method. Table ! provides examples of the generated images.
PLMS4
DDM
05.8
2.00
GHVEL5
GHVB2.5
FLMS4 w HB
GHVB aat
50
71813
GHVB 3.9
L2
8-
0.75
10:
0.50
0.25
80
160
520
640
160
320
640
Number of Steps
Number of Staps Orraal
Figure 10: FID score DIT. Figure 11: 1.2 between different (fours)
Figure 12: Numerical order of convergence for GHVB.
samplers to 1,000-step DDIM.
5.3 EXPERIMENT ON LATENT-BASED DIFFUSION MODELS
We evaluate our techniques using classifier-free guidance diffusion sampling with pre-trained latent- space diffusion models: DiT-XL (Peebles & Xie 2022) and MDT (Gao et (023) (more in Ap- pendix N). We use the -rddr PLMS4 as our baseline for comparison because it demonstrates
superior performance compared other methods (see Appendix L
In Figure 10 significant gap in FID scores on DiT can be observed between 4thh order PLMS4 and 1st ordorder DDIM, this is mostly due to the difference convergence order rather than divergence artifacts. Our GHVB 3.8 and 3.9 techniques successfully mitigate numerical divergence and lead improved FIID scores compared PLMS4, particularly when the number of steps below 10. While HB 0.9 also improves FID scores, using HB 0.8 can worsen FID scores compared to using PLMS4 alone, as shares orrorder convergence with DDIM. For high sampling steps, both HB and GHVB achieve comparable performance to PLMS4 without significant quality degradation. We provide additional results in Appendix L and example images in Figure
5.4 ABLATION STUDY OF GHVB
In this section, we conduct an ablation study of the GHVB method. As explained in Section 4.2 the damping coefficient / of GHVB interpolates between existing AB methods, DDIM and PLMS2. Our goal here to analyze the convergence error of GHVB methods. The comparison is done on Stable Diffusion 1.5 with the target results obtained from 1,000-step DDIM method. We measure the mean L. distance between the sampled results and the target results in the latent space. The results in Figure suggest that the convergence enror of GHVB 1.1 GHVB 1.9 interpolates between the convergence errors of DDIM and PLMS2 accordingly.
Furthermore, we empirically verify that GHVB does achieve a high order of convergence as predicted by Theorem 2 We compute the numerical order of convergence using the formula
log logiewu/ Eau)) where is the error between the sampled and the target latent codes, and & is the q log(kueu/Ke RRaa))
8

---

Published as a conference paper at ICLR 2024
number of sampling steps. As shown in Figure 12 the numerical orders of GHVB 0.5 and GHVB 1.5 approach 0.5 and 1.5, respectively, as the number of steps increases. However, for GHVB 2.5 and GHVB 3.5, the estimated error may be too small when tested with large numbers of steps, and other sources of error may hinder their convergence. Nonetheless, these GHVB methods can achieve high orders of convergence. An analysis of this and other methods is in Appendix
6 RELATED WORK
Accelerating diffusion sampling: Our work is closely related to several other approaches aimed at improving the sampling speed of diffusion models. Certain approaches focus on the improved numerical aspects of samplers, which we have already covered in Section 2.1 Other approaches
involve training additional separate models, including model distillation (Luhman & Luhman [ Salimans & Ho 2022), Schrödinger bridge (De Bortoli et al. 202T), consistency models (Song et al, 2023), rectified flow (Liu et al.] 2022b), and GENIE (Dockhorn et al. (0226). Our goal is an ODE solver that seamlessly integrates with any pre-trained diffusion models. This focus sets our work apart from those aforementioned approaches, which require model retraining.
Our work uses an optimization technique to create an ODE solver. Interestingly, there is rich body of research in the near-opposite direc- tion: using ODEs to analyze optimization methods. Indeed, trajectories made by some optimization methods be described as discrete approximations of ODE solutions. The gradient descent algo- rithm is equivalent to numerically integrating the gradient flow equation: x' = V v=(((())) Polyak (1964) introduced the Heavy Ball method and the Heavy Ball ODE, x" ox' + bV (x())
Connections between momentum methods and ODEs:
O.
which is subsequently investigated by Alvarez (2000); Attouch et al. (2000); Balti & May (016);;
Shi et al. (202T), and [oudoo & Munier (2009). Similarly, Nesterov (1983) introduced the accel-
erated gradient algorithm and derived its associated ODE, given by: x +"++" 7 + V (x((t) O.
Further studies (Su al., 2014; Bubeek et al., 2015 Allen-Zhu & Orecchial 2014 Lessard et al.
2016) analyze similar ODEs to understand Nesterov algorithm better.
Notice that the above momentum associated ODEs are 2nd 2rd-or Nonetheless, we can also regard momentum methods as numerical approximation of the 1 gradient flow, which is very similar to our equation. (Our ODE (2) can be obtained simply by replacing the gradient with the diffusion model's output.) Analyzing our solvers through the -order ODEs is less fruitful as there no guarantee that their trajectories would match those of our DDIM's ODE key goal in our work, inspired by DDIM's excellent quality (De Bortoli et al. 2021
Karras al. 2022).
Momentum in neural network design: The benefits of Heavy Ball ODE, which result in smoother paths and faster convergence to stationary points, have inspired research on neural network archi- tecture design (Moreau & Bruna 2017 Nguyen et al. 2020 Li et al. 2018) and 2odder neural ODE models (Xia et al.] [021]] N Nguyen al. (022). In the context of diffusion models, Wu et al.) (023) reformulate a diffusion ODE into a 2 -order ODE to enhance the training and sampling processes, whereas Dockhorn et al. (2022a) does to the score-based SDE. Instead of redesigning the network with momentum, our paper specifically proposes new momentum basedd solvers for the normal diffusion network to improve its sampling process.
7 CONCLUSION
Our study highlights an issue when using high-order methods diffusion sampling with low sam- pling steps, which leads to the divergence of solutions and artifacts. To address this, we propose two techniques inspired by Polyak's s HB momentum: HB and GHVB. HB serves as an add-on to any diffusion solver to reduce artifacts, while GHVB offers flexibility in selecting the most suitable method for any computation budget. Our extensive experiments confirm the effectiveness of both techniques across image generation tasks with recent diffusion models, such as text-to- text-to-image gen- eration (Stable Diffusion (Rombach et al. 2022)), classifier guidance (ADM (Dhariwal & Nichol 202Th and classifier-free sampling (on DiT (Peebles & Xie 2022), and on MDT (Gao al. (023)

---

Published as a conference paper at ICLR 2024
REFERENCES
Zeyuan Allen-Zhu and Lorenzo Orecchia. Linear coupling: An ultimate unification of gradient and mirror descent. arXiv preprint arXiv: /537, 2014.
Felipe Alvarez. On the minimizing property of second order dissipative system in hilbert spaces. SIAM Journal on Control and Optimization, 38(4): 100-1119, 2000.
Hedy Attouch, Xavier Goudou, and Patrick Redont. The heavy ball with friction method, i. the continuous dynamical system: global exploration of the local minima real-valued function by asymptotic analysis of dissipative dynamical system. Communications in Contemporary Mathematics, 2(01))1-34, 2000.
Mounir Balti and Ramzi May. Asymptotic for the perturbed heavy ball system with vanishing damping term. arXiv preprint arXiv:1609.00135 2016.
Matthew M Berry and Liam M Healy. Implementation of Gauss-Jackson integration for orbit prop- agation. The Journal of the Astronautical Sciences, 2004.
Sébastien Bubeck, Yin Tat Lee, and Mohit Singh. A geometric alternative to nesterov's accelerated gradient descent. arXiv preprint arXiv: 1500xv:500..000 0155506608888 08187, 2015.
Jooyoung Choi, Sungwon Kim, Yonghyun Jeong, Youngjune Gwon, and Sungroh Yoon. ILVR: Conditioning method for denoising diffusion probabilistic models. In 202/ IEEE/CVF Interna- tional Conference on Computer Vision (ICCV), pp. 14347-14356. IEEE, 2021.
Valentin De Bortoli, James Thomton, Jeremy Heng, and Arnaud Doucet. Diffusion Schrödinger bridge with applications to score-based generative modeling. Advances in Neural Information Processing Systems, 34:17695 2021.
Prafulla Dhariwal and Alexander Nichol. Diffusion models beat GANs image synthesis. Ad- vances Neural Information Processing Systems, 34:8780-8794, 2021.
Tim Dockhorn, Arash Vahdat, and Karsten Kreis. Score-based generative modeling with critically- damped langevin diffusion. International Conference on Learning Representations, 2022a.
Tim Dockhorn, Arash Vahdat, and Karsten Kreis. GENIE: Higher-order denoising diffusion solvers. arXiv preprint arXiv:2210.05475 2022b.
Rinon Gal, Yuval Alaluf, Yuval Atzmon, Or Patashnik, Amit H Bermano, Gal Chechik, and Daniel Cohen-Or. An image is worth one word: Personalizing text-to-image generation using textual inversion, arXiv preprint arXiv:2208.0166
2022.
Shanghua Gao, Pan Zhou, Ming- Ming Cheng, and Shuicheng Yan. Masked diffusion transformer is a strong image synthesizer, 2023.
Ian Goodfellow, Jean Pouget-Abadie, Mehdi Mirza, Bing Xu, David Warde-Farles Sherjil Ozair, Aaron Courville, and Yoshua Bengio. Generative adversarial nets. Advances in neural information processing systems, 27, 2014.
Xavier Goudou and Julien Munier. The gradient and heavy ball with friction dynamical systems: the quasiconvex case. Mathematical Programming, 116(1-2)173-1911, 2009.
Desmond J Higham and Lloyd N Trefethen. Stiffness of ODEs. BIT Numerical Mathematics, 33: 285-303, 1993.
Jonathan Ho and Tim Salimans. Classifier-free diffusion guidance. In NewrIPS 2021 Workshop Deep Generative Models and Downstream Applications, 2021.
Jonathan Ho, Ajay Jain, and Pieter Abbeel. Denoising diffusion probabilistic models. In Pro- ceedings of the 34th International Conference on Neural Information Processing Systems, pp- 6840-6851, 2020.
10

---

Published as a conference paper at ICLR 2024
Edward J Hu, Yelong Shen, Phillip Wallis, Zeyuan Allen-Zhu, Yuanzhi Li, Shean Wang, Lu Wang, and Weizhu Chen. Lora: Low-rank adaptation of large language models. arXiv preprint NrXiv:2106009685, 2021.
Xun Huang and Serge Belongie. Arbitrary style transfer real-time with adaptive instance normal- ization. In ICCV, 2017.
Tero Karras, Miik Aittala, Timo Aila, and Samuli Laine. Elucidating the design space of diffusion- based generative models. In NeurlPS 2022 Workshop on Deep Generative Models and Down- stream Applications, 2022.
Tuomas Kynkäänniemi, Tero Karras, Samuli Laine, Jaakko Lehtinen, and Timo Aila. Improved precision and recall metric for assessing generative models. Advances in Neural Information Processing Systems, 32, 2019.
John Denholm Lambert et al. Numerical methods for ordinary differential systems, volume 146. Wiley New York, 1991.
Laurent Lessard, Benjamin Recht, and Andrew Packard. Analysis and design of optimization algo- rithms via integral quadratic constraints. SIAM Journal on Optimization, 26(1):57 95, 2016.
Adam Letts, Chris Scalf, Alex Spirin, Tom Mason, Chris Allen, Max Ingham, Mike Howles, Nate Baer, and David Sin. Disco diffusion. https://github..co com/ altmmicmmmmbc/////mm disco diffusiifi 2021.
Huan Li, Yibo Yang, Dongmin Chen, and Zhouchen Lin. Optimization algorithm inspired deep neural network structure design. In Asian Conference on Machine Learning, pp. 614-629.PMMRR 2018.
Tsung-Y Yi Lin, Michael Maire, Serge Belongie, James Hays, Pietro Perona, Deva Ramanan, Piotr Dollár, and C Lawrence Zitnick. Microsoft COCO: Common objects in context. In Computer Vision-ECCV 2014: 13th European Conference, Zurich, Switzerland, September 6-12, 2014, Proceedings, Part V 13, pp. 740-755. Springer, 2014.
Luping Liu, Yi Ren, Zhijie Lin, and Zhou Zhao. Pseudo numerical methods for diffusion models on manifolds. In International Conference on Learning Representations, 2022a.
Xingchao Liu, Chengyue Gong, and Qiang Liu. Flow straight and fast: Learning generate and transfer data with rectified flow. arXiv preprint arXiv:2209.03000 03300,, 2022b.
Cheng Lu, Yuhao Zhou, Fan Bao, Jianfei Chen, Chongxuan Li, and Jun Zhu. DPM-Solver: A fast ODE solver for diffusion probabilistic model sampling in around 10 steps. In Advances Neural Information Processing Systems, 2022a.
Cheng Lu, Yuhao Zhou, Fan Bao, Jianfei Chen, Chongxuan Li, and Jun Zhu. DPM-Solver++: Fast solver for guided sampling of diffusion probabilistic models. arXiv preprint arXiv22211.01095, 2022b.
James Lucas, Shengyang Sun, Richard Zemel, and Roger Grosse. Aggregated momentum: Stability through passive damping, arXiv preprint arXiv:1804.0033: 2018.
Eric Luhman and Troy Luhman. Knowledge distillation in iterative generative models for improved sampling speed. arXiv preprint axXiv:2101.023380 2021.
Thomas Moreau and Joan Bruna. Understanding the learned iterative soft thresholding algorithm with matrix factorization. arXiv preprint /706.0337 2017.
Yurii Nesterov. A method for unconstrained convex minimization problem with the rate of conver- gence of 12) In Doklady an ussr, volume 269, pp. 543-547, 1983.
Ho Huu Nghia Nguyen, Tan Nguyen, Huyen Vo. Stanley Osher, and Thieu Vo. Improving neural ordinary differential equations with nesterov's accelerated gradient method. Advances in Neural Information Processing Systems, 35:7712-7726. 2022.

---

Published as a conference paper at ICLR 2024
Tan Nguyen, Richard Baraniuk, Andrea Bertozzi, Stanley Osher, and Bao Wang. Momentumran: Integrating momentum into recurrent neural networks. Advances in Neural Information Process- ing Systems, 33:1924-1936, 2020.
Alexander Quinn Nichol and Prafulla Dhariwal. Improved denoising diffusion probabilistic models. In International Conference on Machine Learning, pp. 8162-8171. PMLR, 2021.
Alexander Quinn Nichol, Prafulla Dhariwal, Aditya Ramesh, Pranav Shyam, Pamela Mishkin, Bob Megrew, Ilya Sutskever, and Mark Chen. GLIDE: Towards photorealistic image generation and editing with text-guided diffusion models. In International Conference on Machine Learning, pp. 16784-16804. PMLR, 2022.
William Peebles and Saining Xie. Scalable diffusion models with transformers. arXiv preprint arXiv:2212.09748, 2022.
Boris T Polyak. Some methods of speeding up the convergence of iteration methods. Ussr compu- tational mathematicc and mathematical physics, 4(5):1-17, 1964.
Boris T Polyak. Introduction to optimization. optimization software. Inc., Publications Division, New York, 1:32, 1987.
Alec Radford, Jong Wook Kim, Chris Hallacy, Aditya Ramesh, Gabriel Goh, Sandhini Agarwal, Girish Sastry, Amanda Askell, Pamela Mishkin, Jack Clark, al. Learning transferable visual models from natural language supervision. In International Conference on Machine Learning, pp. 8748-8763. PMLR, 2021.
Saiprasad Ravishankar and Yoram Bresler. Mr image reconstruction from highly undersampled k-space data by dictionary learning. IEEE transactions on medical imaging, 30(5): 1041, 2010.
Daniel Roich, Ron Mokady, Amit H Bermano, and Daniel Cohen-Or. Pivotal tuning for latent-based editing of real images. ACM Transactions on Graphics (TOG), 42(1):1-13, 2022.
Robin Rombach, Andreas Blattmann, Dominik Lorenz, Patrick Esser, and Bjöm Ommer. High- resolution image synthesis with latent diffusion models. Proceedings of the IEEE/CVF Con- ference on Computer Vision and Pattern Recognition, pp. 10684-10695, 2022.
Nataniel Ruiz, Yuanzhen Li, Varun Jampani, Yael Pritch, Michael Rubinstein, and Kfir Aberman. DreamBooth: Fine tuning text-to-image diffusion models for subject-driven generation. arXiv preprint 220iv:2208. 12242, 2022.
Olga Russakovsky, Jia Deng. Hao Su, Jonathan Krause, Sanjeev Satheesh, Sean Ma, Zhiheng Huang, Andrej Karpathy, Aditya Khosla, Michael Bemstein, et al. ImageNet large scale visual recognition challenge. International journal of computer vision, 115(3):211 11553)211-252, 2015.
Tim Salimans and Jonathan Ho. Progressive distillation for fast sampling of diffusion models. In International Conference on Learning Representations, 2022.
Hiroshi Sasaki, Chris G Willcocks, and Toby P Breckon. Unit-DDPM: Unpaired image translation with denoising diffusion probabilistic models. arXiv preprint 2.05552000405358, 2021.
Bin Shi, Simon S Du, Michael Jordan, and Weijie J Su. Understanding the acceleration phc- nomenon via high-resolution differential equations. Mathematical Programming, pp. 1-70, 2021.
Jiaming Song, Chenlin Meng, and Stefano Ermon. Denoising diffusion implicit models. In Interna- tional Conference on Learning Representations, 2020a.
Yang Song and Stefano Ermon. Improved techniques for training score-based generative models. Advances in neural information processing systems, 33: 12248-122448 2020.
Yang Song, Jascha Sohl-Dickstein, Diederik P Kingma, Abhishek Kumar, Stefano Ermon, and Ben Poole. Score-based generative modeling through stochastic differential equations. In Interna- fional Conference on Learning Representations, 2020b.
12

---

Published as a conference paper at ICLR 2024
Yang Song, Prafulla Dhariwal, Mark Chen, and Ilya Sutskever. Consistency models. International Conference on Learning Representations, 2023.
Weijie Su, Stephen Boyd, and Emmanuel Candes. A differential equation for modeling nesterov's accelerated gradient method: theory and insights. Advances in neural information processing systems, 27, 2014.
Xuan Su, Jiaming Song, Chenlin Meng, and Stefano Ermon. Dual diffusion implicit bridges for image-to-image translation. In The Eleventh International Conference on Learning Representa- lions, 2022.
Jinyi Wang, Zhaoyang Lyu, Dahua Lin, Bo Dai, and Hongfei Fu. Guided diffusion model for adversarial purification. arXiv preprint arXiv:2205. /4969, 2022.
Zhou Wang, Alan C Bovik, Hamid R Sheikh, and Eero P Simoncelli. Image quality assessment: from error visibility structural similarity. IEEE transactions on image processing, 13(4);600- 612, 2004.
Daniel Watson, Jonathan Ho, Mohammad Norouzi, and William Chan. Learning to efficiently sam- ple from diffusion probabilistic models. arXiv preprint arXiv:2106.03800 2021.
Daniel Watson, William Chan, Jonathan Ho, and Mohammad Norouzi. Learning fast samplers for diffusion models by differentiating through sample quality. In International Conference Learning Representations, 2022.
Suttisak Wizadwongsa and Supasorn Suwajanakorn. Accelerating guided diffusion sampling with splitting numerical methods. International Conference on Learning Representations, 2023.
Quanlin Wu, Hang Ye, and Yuntian Gu. Guided diffusion model for adversarial purification from random noise. arXiv preprint (tt:::200610877. 2022.
Zike Wu, Pan Zhou, Kenji Kawaguchi, and Hanwang Zhang. Fast diffusion model. arXiv preprint arXiv:2306.0699) 2023.
Hedi Xia, Vai Suliafu, Hangjie Ji, Tan Nguyen, Andrea Bertozzi. Stanley Osher, and Bao Wang. Heavy ball neural ordinary differential equations. Advances Neural Information Processing Systems, 34:1844 1666-6896,,,, 2021.
Qinsheng Zhang and Yongxin Chen. Fast sampling of diffusion models with exponential integrator. In NewrlPS 2022 Workshop on Score- Based Methods, 2022.
Qinsheng Zhang, Molei Tao, and Yongxin Chen. gDDIM: models. International Conference on Learning Representations, 2023.
Generalized denoising diffusion implicit
Richard Zhang, Phillip Isola, Alexei A Efros, Eli Shechtman, and Oliver Wang. The unreasonable effectiveness of deep features as perceptual metric. In Proceedings of the IEEE conference on computer vision and pattern recognition, pp. 586-599 066-595.72 2018.
Wenliang Zhao, Lujia Bai, Yongming Rao, Jie Zhou, and Jiwen Lu. UniPC: A unified predictor- corrector framework for fast sampling of diffusion models. arXiv preprint arXiv:2302.04880) 2023.
13

---

Published as a conference paper at ICLR 2024
Part I Appendices
Contents
Frequently Asked Questions)
15
More Details on Diffusion Modeling and Sampling
16
B.
Diffusion Models and Their Sampling
16
B.2
Numerical Method Summary
17
B.3 Guided Diffusion Sampling
17
Stability Region of Adam- Method
18
The Boundary Locus Technique
19
Derivation of Test Equation
20
Toy ODE Problem
21
Implementation Details of PLMS with HB and GHVB Methods
21
Variance Momentum Methods
21
H Order of Convergence
23
Elaboration on the Order of Convergence Approximation
25
J Additional Qualitative Comparisons on Artifacts Mitigation
25
K Additional Experiment on ADM
26
Additional Experiment on DIT
29
M Additional Experiment on Text-to-Imee Image Generation
29
Experiment on MDT
31
0
Experiment on Personalized Image Generation
32
Face-identity Preserving Image Generation
32
0.2 Reference-Only Control Generation
33
Factors Contributing to Artifact Occurrence in Fine-tuned Diffusion Models
37
Ablation Study HB Momentum
37
14

---

Published as a conference paper at ICLR 2024
R
Ablation Study Nesterov Momentum
42
S
Ablation Study Magnitude Score
42
S.1 Results with Alternative Parameter Settings
43
S.2 Results Alternative Models
43
Statistical Reports
43
A FREQUENTLY ASKED QUESTIONS
Q1: What does the term "divergence artifacts" refer to?
In this paper, the term "divergence artifacts" used to describe visual anomalies that oceur when numerical solutions diverge, resulting in unusually large magnitudes of the results. In the context of latent-based diffusion, we specifically define divergence artifacts as visual artifacts caused by latent codes with magnitudes that exceed the usual range. These artifacts commonly arise when the stability region of the numerical method fails to handle all eigenvalues of the system, leading to a divergent numerical solution. We demonstrate the connection between high latent magnitudes and these artifacts in Figure 13
Scaling factor
Original (x 1.0)
X3.0
x6.0
Stable Diffusion 1.5
Anything Diffusion V4.0
Pastel-mix Diffusion V4.0
Figure 13: Images generated after multiplying the 4 X 4 square at the center of latent codes with constant factor. Using higher scaling factors results more apparent localized divergence artifacts. All samples were generated using PLMS4 (Liu et al.] (022a) with guidance scale of 15 and 250 sampling steps. Prompt: "A beautiful illustration of couple looking at fireworks in a summer festival in Japan".
Q2: How to efficiently determine a suitable solver for diffusion sampling?
Theoretically identifying solver with optimal accuracy and stability region that encompasses all 5A values is challenging because computing À from neural networks is intractable. However, given fixed amount of computational resources (e.g., sampling steps), one gather information about
15

---

Published as a conference paper at ICLR 2024
À to determine the highest solver's order without artifacts by first testing Euler's method with fewer numbers of steps.
In particular, can determine the lowest sampling step of Euler's method that does not produce artifacts (c.g., via a binary search), then use its stability region to determine the highest order of GHVB whose stability region encompasses that of Euler's method. This process does not require additional sampling or evaluating neural networks and can be solved deterministically (Equation
48 5T.. This strategy is more efficient than directly testing different orders the target step because Euler's s method usually requires fewer steps than our target step. For instance, we can infer that GHVB3.2 and those with lower orders are applicable at 20 sampling steps by verifying that Euler's
method works at steps.
Q3: Can we directly interpolate two existing numerical methods instead of using the GHVB method?
Indeed, this possible. However, the order of the resulting method will be the lowest order of the two methods. To illustrate this point, let us consider direct interpolation between the 1 or-order Euler method (AB1) and the 2nd 2rd-order adaass-Bashford method (AB2), expressed as follows:
In+1 In + (1-)))(EEn) (
+ - B - ((((nn)))
(17)
As outlined in Appendix H despite the orders of the interpolated methods, the resulting method is orddrr method.
B MORE DETAILS ON DIFFUSION MODELING AND SAMPLING
This section presents a high-level summary of the theoretical foundation of diffusion models, along with the numerical methods employed in these models. Here, we briefly explain few concepts that contribute to our method.
B.1 DIFFUSION MODELS AND THEIR SAMPLING
Assuming that Fo is a random variable from the data distribution ((zg) that we aim to replicate (e.g., the distribution of natural images), diffusion models define sequence of Gaussian noise degradations of zg random variables Fi I2, .... IT, where
~ q (1) N(
1,
(18)
and Bg e [0, 1] are parameters controlling the noise levels. Utilizing the property of the Gaussian distribution, express 2; directly function of xo and noise NO, I) by
Tt Cr + V1 die,
(19)
where ā, = II 3.). By selecting a sufficiently large T (e.g., T = 1, 000) and an appropriate set of 33, we can assume that 2y follows standard Gaussian distribution (0, I).
The main idea of diffuuson model generation involves training pe (=(x=11/2)) which sequentially re- verses ( using a parametric model while setting p(x;) NO I). Subsequently, we can sample Gaussian noise Fr and use it to evverse-sapll TT- 22, until we retrieve which belongs our data distribution.
The Denoising Diffusion Probabilistic Model (DDPM) Ho et al. 2020) outlines the utilization of neural network denoted as (0 4) to predict the noise e used
computation x, Equation To facilitate the network's training, training image labeled as xp is sampled, alongside values for t and e. These provided values are then utilized in computing 2, using the aforementioned relationship. Consequently, our network is optimized to minimize the difference between the projected noise and the actual noise: ||c- co( t) 12 Once training is completed, we can draw samples _t-1 from the conditional distribution po (a 1x+-1 given by
19
1 - at-11) po po(x =N (1o(5)-1))))
=N
αs 11:
(20)
16

---

Published as a conference paper at ICLR 2024
where μe (z(-11t.t voc (I
VI-a, co( t)) and øT ~ p(ær) = (0, This step-by-step
sampling process is simple but time-consumingg Several sampling techniques and ODE formulations
can be employed to accelerate this process. (See Section 2.1
Many successful diffusion models Dhariwal & Nichol (2021); Song & Ermon 2020); Song et al. (2020b) consider their target distributions go (x as image distributions, with 3y directly repre- senting random variable within the pixel space. These models generate images by inferring net- works directly through diffusion sampling. However, direct training of diffusion models within high-resolution pixel space are computationally expensive. address this, Latent Diffusion Mod- els (LDMs) Rombach et al. (2022) employ a two-stage methodology: Firstly, they employ a trained pair of encoder and decoder D that compress images into smaller spatial latent representations. Secondly, LDMs are then trained on these latent representations
F((=), rather than directly on images æ. Subsequently, the generation of new images involves sampling a latent representation from the diffusion model and decoding it into image using the decoder, D(2).
B.2 NUMERICAL METHOD SUMMARY
Many numerical methods can solve Equations 2 and 3 to accelerate diffusion sampling. In the following section, we will summarize the methods mentioned within our papers. Let's denote our dr f(r), initialized at Tp- ODE as do
Euler's Method is represented by
En+1 = + (xxx),
(21)
where F, represents the approximate solution, and 8 represents the step size. This method, when applied to Equation 2 or results in the DDIM formulation [Snn et al. (2020a).
Heun's Method is 2nd-order extension of Euler's Method given by:
In+1 = In +
2 2 c2):
(22)
Here, = fxx) and e2 = f(x, + 5e1). This method has been used in various diffusion papers, including Algorithm in Karras al. (2022) and DPM- -So--o-- DPM-Solver-2 Lu et al. (2022a), and also serves
as the simplest form of Coorettr methods in Song et al. (2020b): Zhao et al. (2023).
(AB) Methods is family of methods utilizing previous steps to estimate the next step, also known linear multi-step methods. The formulations are follows:
1s- (same as Buler's method):
In+1 = In (f((x))
(23)
nnd orderr
In+1 In + + (3e0 e1),
(24)
-odder:
In+1 = Fn + (23eg - 16e, + 562).
(25)
12
order:
Ent1 = En +
(55cg - 59e; + 37ez - 9eg),
(26)
24
where ek f(xn=R). Notably. SPM- Lu al. 2022b) applies the ndd order AB method
to ODE (3) while PLMS4 Liu et al. (022a)) utilizes the 4 oorder AB method on ODE 2b.
B.3 GUIDED DIFFUSION SAMPLING
Guided diffusion sampling is widely used in conditional sampling, such as text-to- text-to-image and class- to-image generation. There main two approaches for guided sampling:
17

---

Published as a conference paper at ICLR 2024
Classifier guidance (Dhariwal & Nichol
2021
Song al. 2020a) uses pre-trained classifier
model [c Ie, t) to define the conditional noise prediction model at inference time:
E(zr, & | e) (z), sV (c (2,,),
(27)
where 8 "guidance" scale. The model can be extended to accept any guidance function, such CLIP function (Radford et al., (021) for text-lo-image generation (Letts et al.) 2021). This approach only modifies the sampling equation at inference time and thus can be applied to a trained diffusion model without retraining.
Classifier-free guidance (Ho Salimans 2021 trains a conditional noise model cø(2,, è c) to generate data samples with the label c:
E c) = cø(2,, t 6) + s(eg(ze, t c) - ca 6)),
(28)
where & is null label to allow for unconditional sampling. The sampling equations in both ap- proaches can expressed as "guided ODE" of the form
= E(ε, + g(I,o),
(29)
do
where g(z, o) represents a guidance function. To accelerate guided diffusion sampling, splitting numerical methods have been proposed by WWizadwongsa & Suwajanakomi (2023), such Lie- Trotter Splitting (LTSP) divides Equation 29 into two subproblems, i) dy Ely, o) and ii) dz do g(z,ar), but only apply high-order numerical methods to the first equation while resorting to the Euler method for the second equation avoid numerical instability. Higher-order splitting methods, such as Strang Splitting (STSP), can also mitigate artifacts. However, these methods require solving the second equation twice per step, which comparable to increasing the total sampling step to avoid artifacts. Both approaches require non-negligible computation.
dg
C STABILITY REGION OF
METHOD
To investigate the stability of the AB2 method, we apply AB2 to the test equation x' = Ax, which 1 11nn-1) To 2
was also used with the Euler method (Section 2.2p. We have 2n+1 = In + (āAr, - solve this linear recurrence relation, we substitute = r" into the formula, where is complex constant. Simplifying the resulting equation, we obtain the characteristic equation:
3 θf
1 06 0,
2
1
(30)
2
2
which has the solutions
r1
1 1+
3 2
82
v(1+
V
|+
3 2
5d
28A
(31)
r2
1+ 82
3
v(1
V
1
+
3
3
82
28A
(32)
2
2
The general formulation of n, can be expressed as
En airl + 027
(33)
where aj and d2 are constants. The numerical solution 2n, tends to 0 as n tends to infinity when both Iril < 1 and |r2l < 1, which means the stability region of AB2 is determined by the complex region
((1
(1+
s = z EEC:
+ 2
=
2z
(34)
2
Solving for the complex area from the roots of the characteristic equation can pose significant chal- lenges in numerical analysis. One commonly employed graphical technique to visualize the stability region is the boundary locus technique (Lambert et al. 1991).
18

---

Published as a conference paper at ICLR 2024
C.1 THE BOUNDARY LOCUS TECHNIQUE
The boundary locus technique begins by defining the shift operator E such that Exg = -.- Note that EaIa= Eag-1 2 Generally, numerical method can be represented in the following form:
a (E.))
(35)
where A and are polynomials of E. For example, in the case of the AB2 method, we have A(E) = - E and B(E) -
To determine the stability region of numerical method, we apply the boundary locus technique to the general form given by Equation 35] The characteristic equation of the method can be obtained by substituting ((x)) Arn (i.e., the test equation) and iin = r", which yields
A(r- ())(r-)),
(36)
where the root of the method's characteristic equation. The stability region of the method is the area in the complex plane where the characteristic root have modulus less than 1. The
boundary of the stability region can be determined by substituting r with a modulus of (which means that e for some real value 9) into the characteristic equation and solving for 8X. This yields the locus of points in the complex plane where the characteristic roots of the method the boundary of the stability region. Specifically, B(e-(e) A(e=t))/B(ee
can obtain the curve = s(0) where E [-x, v). that represents the boundary of the stability region in the complex plane. By comparing the stability regions of different numerical methods, we can determine which method is more stable and accurate for a given problem. The boundary locus technique provides powerful tool for analyzing the stability of numerical methods and can help guide the selection of appropriate methods for solving ODE problems.
A(e-1)/
Example C.1. (Euler Method) The Euler method, a numerical technique for approximating solu- tions of ODE, can be expressed as:
E)= = 8Ef(z))
(37)
The associated polynomials for this method are:
A(=) 1- B(2)= ==
(38)
The stability region of the Euler method corresponds to the locus curve in which the solution remains bounded. This region can be determined by evaluating the complex function:
Ale-i0
e-e-
s (0)
1-1,T
(39)
B(c-i8)
c-i6
The locus curve forms perfect circle with a radius of 1 and a center at -1. Example C.2. (AB Methods) The 2nd- -order Adams- Adams-Bashpee Bashforth (AB2) method is given by:
(1- -)E, = E- (2
f(xx).
(40)
The locus curve representing the stability region of this method is given by:
1 e
2(1 -00
0 E 1-
(41)
250
3e 18
226
2
Similarly, the stability regions for the AB3 and AB4 methods can be obtained by evaluating the complex functions:
12(1 e -10
s(0)
1 O [E[ T..
(42)
23e-id i0
16e 228
+ 55---00 3i0
24(1 - -i0)
s(0)
0 E [-
(43)
55c
59c + 37e 9E
The locus curves for the boundary of stability regions of the first four AB methods are visualized in Figure
19

---

Published as a conference paper at ICLR 2024
Additionally, include the locus curves for the boundaries of the stability regions of our tech- niques. As shown in Figure the local curves for the PLMS method with HB B are given by: PLMS1 with HB 8:
s(0) (1 e e- (1 B)e-10)
(44)
9e-10 10
PLMS2 with HB 3:
s0)) 2(1 - e e010) (1 - (1 - 8) 90)-100)
(45)
8(3e e- 210) e2200)
PLMS3 with HB 8::
s(0) = 12(1 e 18 (1 (1 3) 18)
(46)
B(23c- 18 16c + 53-310
PLMS4 with HB 8::
24(1 e -00 (1 B) :0)
s0) 355e 10 59e-59--0 210
(47)
+ 37e
310 9e- 1:0)
Similarly, the local curves for the GHVB method in Figure 6 are given by: 1 GHVB (equivalent to PLMS1 with HB):
(1 E (1 (1 - 3)e
8 (0)
(48)
Be- 18
2nd -oddr GHVB:
s(0) 2(1 - c e010) - (1 - 3)e
(49)
((2 8)e (2 3)e 210)
3rd -order GHVB:
12(1 - e 10 (1 (1 B)E -00)
s(0)
(50)
(18 + 55)e 10 (24 8.3e 260 (6 B)e 330
4th -rder GHVB:
1 (1
B) =03
s(0) (46 + 9:3)--0- (78 - 193)e-216 210 + (42 5,3 e 320 (10 - B)e 10ee-100 128
(51)
D DERIVATION OF TEST EQUATION
This section presents the derivation of the test equation v' Au,, which serves as a fundamental tool for analyzing the stability of numerical methods in diffusion sampling, discussed in Section 3.2
Starting from the differential equation for F. we have
dr do Vé(x")(=-1))
(52)
We then define u vJ - ))), where u a normalized eigenvector of Vē(a*)? corresponding to the cigenvalue X. Taking the derivative of with respect to a and using the chain rule, we have: du do do (o -'==')=='iV((p> (V(("""")) x²)
T d
(53)
(u) (Av) x²)
(54)
Au
(55)
Thus, obtain the test equation u' = Au.
20

---

Published as a conference paper at ICLR 2024
E TOY ODE PROBLEM
This section aims to demonstrate how the solutions yielded by numerical methods can diverge when the stability regions of the methods too small. Additionally, we illustrate how momentum- based techniques can enlarge the stability region. The demonstration is conducted on a 2D toy ODE problem given by:
dx
0
1
11
-
(56)
dt
10
02
The eigenvalues of the 2 x 2 matrix are -9 and -1, and the exact solution of Equation 56 is given by:
1
<t))
98
(57)
8
8
As increases, rtf) converges to the origin.
Let us say we want to numerically compute x²3) by integrating the ODE for 26 steps with nu- merical method. We divide the time interval [0, 3] into 26 equal intervals, resulting a step size of 5 0=33 3/26. For this particular setting, it turns out that the 2nd diverges, but the Fuler method converges. To see this, observe that the stability region of the AB2 method only cover the interval |-1, 0] of the real line, as depicted in Figures [4b] and [4d So, for the eigenvalue À -,, the product dd = 227/26 lies just outside the region. Consequently. the numerical solution yielded by AB2 diverges, indicated by the blue line in Figures 14a and 14c In contrast, the Euler method's 's stability region contains both values of 5λ, and the numerical solution, represented by the green line in Figures 14a and [4c appears be more accurate.
Adams-Bashforte (AB2) method
We can improve the stability from AB2 method by applying any of our proposed techniques: Heavy Ball momentum (HB) and GHVB. The stability regions of the modified AB2 methods given by the red (3 0.8) and yellow (3 = 0.9) lines in Figures [4b] and [4d Observe that they contains the points associated with the AO values. As result, the numerical solutions of the modified
methods converge to the origin, as demonstrated by the red and yellow lines in Figures [4a and 14c respectively.
(a) solution
(b) stability region
(c) solution
(d) stability region
Figure 14: Comparison of solution trajectories and stability regions of various numerical methods when applied to the toy ODE problem. Here, we seek compute xx3) in 26 steps with the Euler method, the AB2 method, and methods resulting from modifying AB2 with our momentum-basee techniques. Subfigure (a) presents the numerical solutions obtained using our modified AB2 method with HB momentum, while subfigure (c) showcases those obtained using our GHVB. The stability regions of the methods are depicted in subfigures (b) and (d) respectively.
F IMPLEMENTATION DETAILS OF PLMS WITH HB AND GHVB METHODS
In this section, we present the complete algorithms for the PLMS method with the HB momentum and the GHVB method in Algorithm and Algorithm [ respectively.
G VARIANCE MOMENTUM METHODS
In 2018, a variant of Polyak's HB momentum called aggregated momentum (ucas et al. 2018) was proposed. Its objective is enhance stability while also offering convergence advantages.
21

---

Published as a conference paper at ICLR 2024
Algorithm 1: PLMS step with HB momentum
input: 7. (previous result), (step size), (er (evaluation buffer), r (method order),
Pn (previous ;
en = Eo (Tx):
= min(r,n)
if c = then
= Cn
else if = 2 then
= (3ea en- /
else if c = 3 then
è (23en 16en-1 +
else
é = (55en -1+ 3Ten--3 - Pen 2)244
Un+1 (1 B) + -Ø;
Result: In + dun+1
Algorithm 2: GHVB step
input: nn (previous result), (step size), 3 (damping parameter)
(usiign (evaluation buffer), r (method order),
on+1 (1 B)E, + (Z,))
mi)(r,n)
if c ------ then
= Pa+1:
else if == 2 then
= ((2 + 3) (2 B)un)/23
else if c = 3 then
è= ((18 + 6,3) - (24 88)
+(6 8) 1)/1///1122
else if c = 4 then
((46 9.3) 9.9)) un+1
(78 193)un
+(42 5,3) 11
(10 ())nn=))220 2/228
else
(3420- 646,3) Un + (2880- 264 Un (1380- 106,3)nn22
((1650 + 251,3)nn.)
+(700 193)u,,a)
Result: I, + dē
This modification introduces multiple velocities, denoted by (6) each associated with its specific damping coefficient g())
K
n (1 - ((()) + ((xn), In+1 = Iin +
(58)
Nesterov's momentum (Nesterov) 1983 is version of the classic momentum that can also be applicd to diffusion sampling processes to improve stability. It can be written as follows:
yn+1 = In + (zx()))
In+1 Dn+1 + (1 B) -y))
(59)
In fact, Nesterov's momentum can be obtained from aggregated momentum by considering the following:
# (1 - ())))) + +f((x,), Jn+1 = 2n + 5((1 - -3) (22))
(2) = f(xn),
n+
(60)
The stability regions of Nesterov's momentum when applied to the Euler method and high-order tdams-Bashfore methods illustrated in Figures 15a through 15d Observe that the stability
22

---

Published as a conference paper at ICLR 2024
regions of methods with Nesterov's momentum become larger in similar manner to those with Polyak's 's HB momentum. However, the enlargement due to Nesterov's 's momentum is more pro- nounced in the vertical direction, while the Polyak's HB momentum's 's enlargement is more horizon- tal nature. (See Figures and 14d) The differences the shapes of the stability regions suggest type of momentum is more suitable certain ODE problems than the other.
10-
60
-65
2010
433
-03
-15
-25
-25 -22
(a) 1d order
(b) 2nd order
(c) 3 order
(d) 4th order
Figure 15: Comparison of stability regions for different methods with different levels of Nesterov's momentum.
While generalizing the aggregated momentum method to higher-order methods is possible, is no longer straightforward as it is with the HB method. As an example, we will consider the 2rd-order generalization of Nesterov's 's momentum method.
We begin by noting that (3 (1 3)A)un+1 = (f(xx)) Our goal is to find polynomials B and C such that
Ar 8(B(A)en+t + C(A)f(zn)).
(61)
Multiplying both sides of the equation by (3 (1 3)A). get
(8 + (1 B) ((B (z,) (1 (1)))))(()))))).))))
(62)
Replace the left side with the first two terms from Equation 13
we obtain
8 (B+ 2 2 )) f(x,) = 8(BB(A) ((BB(A)(BBB)))
(63)
2
Let B(A) bn b,A and ((A) = Co. Then, by balancing the coefficients A on both sides, we
have 1 = bo + Co and 2+8
2 = Bb, + (1 -B) We can now write the final formulation as follows:
In+1 on + 8(b08n+1 + bj (fn+1 - vn) + cof))
(64)
This suggests that there are countless different ways to expand the stability region of a numerical method, which offer many new research opportunities.
H ORDER OF CONVERGENCE
When solving ODEs numerically, it is important to consider the accuracy of the method used. One way to measure accuracy is by considering the order method's of convergence the. Suppose we have numerical method of the form
A(E),, (E())((E))
(65)
where A(E)- + a + un E2 +... a E and B(E) by + b, + E², The method is said to be of th order if and only if, for all sufficiently smooth functions x, we have that
an - mō) -
bm (a m5) = (
(66)
m=0
where x' denotes the derivative of s..
23

---

Published as a conference paper at ICLR 2024
To derive the order of convergence,
Taylor expansion for both and x', yielding
md)
md)*
L.H.S.
Um
(1)
(0)
5
bm
(o)
k!
k!
m
k=0
m
k=
(-mo)
x (k) (0)
Dm
x x((+1)(0)
k!
k!
k=0
m==
k==
m==
(-m8)h
m m ( (-0)
x.
(8) (a)
bm
x²
k!
(h: 1)!
k=
mm=0
R=1
mm=
m
x (o)(
mm
(k:
1}!
k=
where (a) denotes the htth derivative of x evaluated at a. Therefore, the method has ndh order of convergence if and only if the coefficients satisfy the condi- tions given by
am
m
772
k-1
bm
(k 1)! 0,
0,
(67)
k!
Now, we discuss the convergence order of any numerical method after HB momentum is applied to it. An example of such an algorithm is the modified PLMS method, presented Algorithm
Theorem (Convergence order of numerical methods with HB momentum). Suppose that a pmh
order mumerical method has the form In+1 In + 8 5 bm If (In-m). where 2 1. The modified method that uses HB momentum can be expressed as follows:
Cn+1 (1 -1- -gun +
Σ (#n-m)
(68)
m==
In+1 IIn + 8un+1.
(69)
has first-order convergence.
Proof. From the condition given in 67 it follows that Em by = We can rewrite these equa- tions as:
In+1 Zn. (1 (2,, xn-1) 83 bm (I
(70)
m1=00
To estimate the order of the modified method, we evaluate Equation 70 and obtain:
am = - = 0;
(71)
m==
am
m'
bm m
0 (1 () 2)
Σ bm
0!
= 3 0.
(72)
Therefore, we have shown that the modification to the method has first-order convergence.
L
Next, we turn our attention to the GHVB method.
Theorem 2 (Convergence order of the GHVB method), The th order GHVB (Algorithm order of convergence of
has
24

---

Published as a conference paper at ICLR 2024
Proof. We will use the ordde method as example. Using Equation 16 can write an equiv-
alent equation as:
2n+1 Zn -1 -1-1111111 (I xn-1) = 2 2 f(zn)
2-3 3 f(aa-i))
(73)
2
To estimate the order of the modified method, we evaluate Equation 70 and obtain:
and
Om
1 - (1 - B1 - =0
m=0
m
(1- 2) 22+
2 8
Om
bm
0
0,
1!
0!
2
Qm
m
2
Σ
bm
m
1 (0
00(11-290 - 22 1 (1- -11
(2+ B
1
2 - 3 2
2!
1!
2
2
m=0
m==0
2
3,3
2
33
O.
2
Thus, the method has convergence order of two. Methods of other orders can be dealt with in similar fashion.
I ELABORATION ON THE ORDER OF CONVERGENCE APPROXIMATION
In Appendix H. we explored the theoretical aspects of the order of convergence for numerical meth- ods. In this section, we will delve into the estimation of the order of convergence specifically for GHVB Section 4.2
To assess the order of convergence, we focus the error e, referred to the global truncation error. This error quantified by measuring the absolute difference the latent space between the numerical solution and an approximate exact solution obtained through 1,000-step DDIM sampling. The order of convergence for a numerical method is defined as g, where the error e follows the relationship (8)), with å representing the step size.
To estimate the order of convergence practically, we adopt a straightforward approach. It involves selecting two distinet step sizes, denoted as &new and dold. and computing the corresponding emross Enew and Gold- These errors can be approximated using the following formulas:
Cnew Rd Cnew (Onew) )²;
Cold R Cota(oata)
(74)
Here, we make the assumption that Cnew is approximately equal to Cola- By taking the ratio of Cnew to Cold. we obtain:
Cnew
8pew
(75)
Cold
3old
Consequently, we can estimate the order of convergence, denoted as q, by evaluating the logarithmic ratio of errors and step sizes:
log(Enew/ / (Cota)
q log( Loge / (Saa)
(76)
In our investigation of GHVB in Section 4.2 we conducted sampling experiments using 20, 40, 80, 160, 320, and 640 steps. This choice of an exponential sequence for the number of steps was inten- tional, it allowed to approximate / / 1/2. By doing so, facilitated the estimation process. The results, representing the approximated order of convergence for GHVB, are visually depicted in Figure 12
ADDITIONAL QUALITATIVE COMPARISONS ON ARTIFACTS MITIGATION
Figure compares momentum- based methods, HB and GHVB, with two different diffusion solver methods, DPM-Solver++ (Lu et al.) 20226) and PLMS4 (Liu et al., [022a), without mo- mentum. The number of sampling steps is held constant while varying the guidance seale to
25

---

Published as a conference paper at ICLR 2024
intentionally induce divergence artifacts. (Note that the guidance scales that yield such artifacts different between diffusion models.) The figure demonstrates that, under the difficult settings low step counts and high guidance scales where the baseline methods produce artifacts, proposed techniques can successfully eliminate them.
We present additional qualitative results to show the effect of the damping parameter the quality of images generated by methods modified with HB momentum. We use methods of varying orders, including DPLM-Solver++ (Lu al., (022b). UniPC (Zhao et al. 2023), and PLMS4(Liu al.)
2022a). The diffusion models utilized in our analysis are Realistic Vision v2.0 Anything Diffusion
v4.0 Counterfeit Diffusion V2.5 Pastel-Mix Deliberate Diffusson and Dreamlink Diffusion
V1.0%
The results are shown in Figure 16 Notice that stronger momentum (lower 3) leads to fewer
and less severe artifacts.
K ADDITIONAL EXPERIMENT ON ADM
We present additional details and results for the ADM experiment Section 5.2 The primary objective of this experiment was to provide quantitative evaluation of class-conditional diffusion sampling in the context of pixel-based images. The experiment was conducted using the pre-traines diffusion and classifier model at the following link: 7 The implementation used in our experiment was obtained directly from the official DPM-Solver GitHub repository
To enhance the capabilities of DPM-Solver++ we incorporated HB momentum into DPM-Solver++ (just change few lines of code) and implemented the splitting method LTSP both with and without HB momentum into the DPM-Solver code for comparative purposes. The experiment was done four NVIDIA RTX A4000 GPUs and a 24-core AMD Threadripper 3960x CPU.
The results of the experiment, measured by the full FID score, are presented in Table (as well as in Figure 99. Our technique is highlighted in grey within the table, while the best FID score for each number of steps indicated in bold. Additionally, Figure 17 showcases sample images from this experiment. As this experiment uses pixel-based diffusion models, the observed divergence artifacts differ from those latent-based diffusion models. Specifically, the artifacts may display excessive brightness or darkness caused by pixel values nearing the maximum or minimum thresholds.
Number of Steps
Method
12
14
16
18
20
25
DPM--Solve+++
66.77
46.77
34.56
26.97
21.87
19.48
16.31
15.63
DPM-Solver++ w/ HB 0.8
47.10
33.65
25.61
21,42
18.94
17.76
15.98
15.53
LTSP [PLMS4, PLMS1]
45.32
34.08
26.58
21.54
18.54
17.15
15.79
15.51
LTSP [PLMS4, GHVB 0.8]
37.43
29.74
23.60
20.23
18.46
17.14
16.07
15.79
Table 1: FID scores on classífier-guidanee ADM models
httes///hnnginneeeec/////22///21 Vision V2.0
https: httear//angiigaeeee/////////oee2
26

---

Published as a conference paper at ICLR 2024
Without Momentum
With IIB 0.8 (Ours)
With HR O.5 (Ours)
Realistic Vision V2.0
Anything (4.0 20)
Deliberate is
(a) DPM (2lMer++(M))) using 15 steps and guidance scale Prompt: "a tiny cute bunny"
Withous
With (Ours
With (Oum)
Anything (4d
Pastel-nix V4,0(e
Dicamlike (1.0
(b) UniPC using 8 steps and guidance scale Prompt: "a kitty"
Momentum
(Ous)
(Our)
V1.0
(c) PLMS4 using steps and guidance scale s. Prompt: "cute humanoid red panda"
Figure 16: Impact of different damping coefficients 3 on HB momentum for 2 DPM- Solver++(2M) (Lu et al. 2022b).. 3rd-order UniPC (Zhao et al. 2023). and 4th-order (Liu et al., 2022ap. Incorporating higher momentum values (lower (3 helps mitigate the occurrence of diver- gence artifacts.
27

---

Published as a conference paper at ICLR 2024
DPM-Solver++
DPM-Solver++ with HB 0.8 (Ours)
LTSP [PLMS4, PLMSI]
LTSP [PLMS4, GHVB 0.8] (Ours)
Figure 17: Samples from DPM-Solver++ (Lu et al. 2022h and LTSP (Wizadwongsa & Suwa- janakorn 2023) with and without momentum, employing classifier guidance diffusion with a guid- ance scale of 10 and 20 sampling steps.
28

---

Published as a conference paper at ICLR 2024
L ADDITIONAL EXPERIMENT ON DIT
This section supplements the experiment on DiT in Section 5.3 and reports results for additional baselines; DPM- + (Lu et al., 2022b) and LTSP4 Wizadwongsa & Suwajanakorn 2023 The code implementation and pre-trained DiT model were obtained directly from the official GitHub repository g
The experiment was done on four NVIDIA GeForce RTX 2080 Ti GPUs and a 24-core AMD Threadripper 3960x CPU.
Table presents the full results including those of the original baselines. Our technique is high- lighted in grey in the table, and the best FID score for each number of steps is in bold. We also include the improved Precision and Recall metrics (Kynkäänniemi et al. 2019) in Tables 3 respectively, where higher values indicate superior performance. Generated samples for different sampling methods are shown in Figure [8] Our methods produce comparable FID scores LTSP4 and PL.S4 at 25 steps, but outperform all baselines for lower 6-20 steps.
Number of Steps
Method
10
15
20
25
DDIM
55.35
36.97
26.06
19.47
15.02
8.04
6.52
5.94
DPM-Solver++
18.60
10.80
7.93
6.72
6.13
5.49
5.30
5.24
LTSP4 [PLMS4, PLMSI]
13.33
9.01
7.49
6.55
6.09
5.32
5.20
5.17
PLMS4
13.10
8.94
7.31
6.51
6.03
5.32
5.21
5.17
PLMS4 w/ HB 0.8
14.35
9.25
7.46
6.68
6.19
5.47
5.29
5,24
PLMS4 w/ HB 0.9
11.66
8.16
6.69
6.21
5.78
5.29
5.19
5.17
GHVB 3.8
10.99
7.93
6.63
6.19
5.80
5.31
5.22
5.18
GHVB 3.9
11.67
8.29
6.83
6.30
5.87
5.31
5.22
5.18
Table 2: FID scores on DiT-XL
Number of Steps
Number of Steps
Method
8
10
20
Method
6
10
20
DDIM
0.36 0.56 0.67
0.79
DDIM
0.60 0.64 0.65
0.67
SPM-Solver++
0.63
0.75
0.79
0.81
DPM-Solver++
0.67 0.68
0.68
0.68
LTSP4
0.67 0.74 0.78
0.81
LTSP4
0.70 0.70 0.69
0.68
PLMS4
0.68 0.75 0.78
0.81
PLMS4
0.70 0.70 0.69
0.68
PLMS4 w/ HB 0.8
0.68 0.77 0.79
0.81
PLMS4 w/ HB 0.8
0.68 0.68 0.67
0.68
PLMS4 w/ HB 0,9
0.70 0.77 0.79
0.81
PLMS4 w/ HB 0.9
0.69 0.69 0.67
0.68
GHVB3.8
0.71 0.77 0.79
0.81
GHVB3.8
0.69 0.70 0.68
0.68
GHVB3.9
0.70 0.76 0.78
0.81
CHVB3.9
0.70 0.70 0.69 0.68
Table 3: Precision DiT-XL
Table 4: Recall on DiT-XL
M ADDITIONAL EXPERIMENT ON TEXT-TO-IMAGE GENERATION
To provide more comprehensive evaluation of the methods discussed in Section 5.1
we utilize fine-tuned variant of Stable-Diffluion called Anything V4. We consider full-path samples generated by PLMS4 (Liu al.) 2022a) at 1,000 steps reference solutions. The performance of each method evaluated by measuring the image similarity between the generated samples produced using a reduced number of steps and the reference samples. Importantly, both sets of samples originate from identical initial noise maps. This comparison allows us to assess how well the solution from each configuration matches the full-path reference solution.
29

---

Published as a conference paper at ICLR 2024
Number steps
25
PLMS4 Ciu ot
PLMSA
GHVR338
Liee etal PUZ22
PLMS4 w/ TIB 0.9
GHIVB388
Comparison of samples generated from DiT-XL with guidance scale of 3, using differ- ent sampling methods and sampling steps.
Figure 18:
We measure image similarity using Learned Perceptual Image Patch Similarity (LPIPS) (Zhang jet al. (0188) and the L2 norm in the latent space, as discussed in Section 4.2 The outcomes of these analyses visually presented in Figure [9] We also include the results from other metrics, e.g.. Structural Similarity Index (SSIM) (Wang et al. 2004) and High-Frequency Error Norm (HFEN)
(Ravishankar Bresler 2010). in Figure 20
0.09
1.30
023
64-
0.33
Figure 19: Comparison of LPIPS and L2 distance across different sampling methods, with and without using our momentum techniques. The experimental setting is similar that in Section 5.1
30

---

Published as a conference paper at ICLR 2024
056
GIN022
OPARO
50t
0.6-
S1FG2
010R0
GHA2L
0.4
Number Staps
Number Slags
Figure 20: Comparison of SSIM (higher is better) and HFNE (lower is better) distance across differ- ent sampling methods, with and without using our momentum techniques. The experimental setting similar that in Figure II
Discrepancies between the numerical solutions and the 1,000-step reference solution arise from two primary factors: the accuracy of the employed method and the presence of divergence arti- facts. Notably, in this particular context, divergence artifacts tend outweigh errors stemming from method accuracy. Consequently, higher-ordet methods exhibit greater discrepancies in both LPIPS and L2 similarity measurements. able ability to minimize deviations from the reference solution in the majority of cases. Furthermore, Figure 19 consistently demonstrates the superiority of our techniques compared to other methods, as evidenced by the lower LPIPS and L2 similarity scores. These results indicate that our techniques effectively both reduce divergence artifacts and handle errors related method accuracy.
worth highlighting that our techniques demonstrate remark-
N EXPERIMENT ON MDT
In this section, conduct additional experiment using Masked Diffusion Transformer (MDT)
Gao et al. (2023) for classifier-free guidance diffusion sampling. The code implementation and
MDT model pre-trained on ImageNet256 (Russakovsky et al. 2015) were obtained directly from
the official GitHub repository [10 The experiment was done on two NVIDIA GeForce RTX 3090
GPUs. Similar the experiment on DīT-XL (Peebles & Xie 2022) in Section [.3] we compare FID
scores of numerical solvers with our proposed momentum techniques to the baseline of PLMS4 (Liu
al. 2022a) and DDIM (Zhang et al.) 2023) using 50,000 generated samples.
As presented in Table both HB and GHVB momentum can mitigate artifacts and lead to improved FID scores compared PLMS4. Note that using HB 0.8 yields worse FIID scores compared using
PLMS4 alone because of its 1 order convergence property. Similar to Appendix L we also include
the improved Precision and Recall metrics (Kynkäänniemi et al.) 2019) in Tables X and 7
This
experiment further highlights that our methods work across various recent diffusion models.
Number of Steps
Method
6
7
8
9
10
15
20
25
DDIM
56.53
36.79
24.17
17.13
12.62
5.11
3.40
2.77
PLMS4
11.10
6.44
5.45
4.30
3.89
2.53
2.15
2.01
PLMS4 w/ HB 0.8
13.36
7.19
5.03
4.03
3.50
2.42
2.12
2.03
PLMS4 w/ HB 0.9
10.48
5.77
4.38
3.65
3.26
2.28
2.03
1.96
GHVB 3.8
9.67
5.60
4.48
3.80
3.39
2.36
2.07
1.98
GHVB 3.9
10.18
5.88
4.83
3.98
3.58
2.43
2.10
1.99
Table 5: FID scores on MDT
"https::/githubb.cm 1
31

---

Published as a conference paper at ICLR 2024
Number of Steps
Number of Steps
Method
6
8
10
20
Method
6
10
20
DDIM
0.33
0.53
0.64
0.76
DDIM
0.50
0.54
0.57
0.63
PLMS4
0.66
0.71
0.73
0.76
PLMS4
0.62 0.65
0.66
0.66
PLMS4 HB 0.8
0.64
0.72 0.74 0.77
PLMS4 w/ HB 0.8
0.59
0.63
0.64
0.64
PLMS4 w/ HB 0.9
0.66
0.73 0.74 0.77
PLMS4 w/ HB 0.9
0.61 0.64
0.65
0.65
GHVB3.8
0.67 0.72 0.74
0.77
GHVB3.8
0.62
0.65
0.66 0.65
GHVB3.9
0.67 0.71 0.73 0.77
GHVB3.9
0.62 0.65 0.66 0.66
Table 6: Precision on MDT
Table 7: Recall on MDT
O EXPERIMENT ON PERSONALIZED IMAGE GENERATION
This section evaluates our methods on personalized image generation, a task which involves gen- erating images with specific concepts, styles, or following certain compositions given by reference images. The sampling process in this task is significantly more challenging to speed up due to greater incidence of artifacts compared to the regular diffusion sampling conditioned with only a text prompt. All experiments were conducted on four NVIDIA RTX A4000 GPUs.
0.1 FACE-IDENTITY PRESERVING IMAGE GENERATION
This experiment focuses on generating images of people's 's faces with certain face orientations, speci- fied by face key points. We fine-tuned Stable Diffusion (Rombach et al. 2022) using 6 rank-1 LoRAs (Hu et al. 2021 across six facial identitics, which was implemented in the GitHub repository We
|1|
used face dataset accessible at Kaggle For each identity in the dataset, the corresponding LoRA
was trained with 10 randomly selected images using Pivotal Tuning (Roich et al.] 2022). ) consisting of Textual Inversion (Gal et al. 2022) followed by LoRA DreamBooth (Ruiz al. 2022) training,
for 1,000 steps. Note that we employed Stable Diffusion 1.5 L33 as the base model. Our proposed HB methods were applied to PLMS4 (Liu et al. 2022a) and D (Lu al. 2022b), which also served as baselines for our comparison.
At test time, we utilized a variant of Stable Diffusion named DreamShape and applied the LoRAs with fixed LoRA scaling parameters a = 1.0 and a = 0.4 for the text encoder and the UNet model, respectively, resulting in different models. The sampling process is further conditioned with the outputs from the ControlNet V1. V [rrinndd on a human pose estimation task, Subsequently, we prepared 25 sets of inputs, i.e., text prompt, face key points, random seed, and a fixed guidance scale of 15, and generated 25 converged, reference results by applying 1,000-step PLMS4 sampling to each of the augmented models. resulting in a total of 150 target images. For evaluation, we calculated the mean LPIPS ZZhng al.) 2018) and L2 distance between these images latent and those generated by each ODE solver when provided with the same inputs.
For each numerical solver employing momentum number 3, we conducted grid search find the optimal value in the range (0.1,...,0.9) 0.9] that yields the least number of sampling steps such that the LPIPS comparable that of 200-step DDIM Zhang et al. (2023). As shown Figure 21] the results indicate that, with the same number of sampling steps, PLMS4 with HB and GHVB generated higher-quality images than vanilla PLMS4. This implies that those with momentum can achieve the target image quality with fewer steps. In particular, our two methods reached the LPIPS threshold set by the 200-step DDIM method with 16.53% and 25.38% less sampling time, respectively, compared to PLMS4 without momentum. A similar result can be observed when comparing +PM-Solver++ and its HB momentum variant, resulting in 18.22% less sampling time. Figures 23 and 24 provide examples of the generated images.
httpp:://huppn
5_opeppoee
32

---

Published as a conference paper at ICLR 2024
PLMS4 HB 0.8
PIMSA wy HB D.t GHVBS.1
GRVB3.1
DPN-Solowre+ HB
Moun La Distancee
DPM-Scleere-
UPD5
2x10H1
150
50
100
(a) LPIPS
(h) Mcan 1.2 distance
Figure 21: Comparison of LPIPS and 1.2 distance of (a) PLMS4 with and without momentum, (b) DPM-Solver++ with and without momentum, and (c) GHVB in the face-identity preserving image generation experiment.
O.2 REFERENCE-ONEE CONTROL GENERATION
We evaluate our techniques on the reference-one control pipeline implemented for Stable Diffusion in the following GitHub repository 16
This pipeline conditions the sampling process text prompt and reference image, enabling content and style transfer of the reference image via attention mech- anism and feature rescaling proposed in AdaIN (Huang & Belongie 2017
We utilize a fine-tuned Stable Diffusion model named Anything V3.0 to generate Japanese animation-style images, Specifically, prepared 150 sets of inputs, i.c., text prompt, reference image, random seed, and fixed guidance scale of 7.5 and style fidelity of 1.0 (see ( and gener- ated the target results by using 1,000-step PLMS4 (Liu et al. 2022a). For evaluation, we calculated the mean LPIPS (Zhang et al. 2018) and L2 distance between the couve-converged result and those generated by numerical solvers when provided with the same inputs.
The results in Figure 22 suggest that PLMS4 with HB and GHVB generated higher-quality images than the regular PLMS4. Specifically, our two methods reached the LPIPS threshold set by the 150- step DDIM method with 22.40% and 42.92% less sampling time, respectively, compared to PLMS4 without momentum. Figure 25 provides examples of the generated images.
L2 Distance
3x12
20*
LPPS
2x10-1
10N
100
150
200
250
300
50
100
155
200
250
300
Number staps
Nambbet stags
(a) LPIPS
(b) Mean L.2 distance
Figure 22: Comparison of LPIPS and L2 distance of (a) PLMS4 with and without momentum and (b) GHVB in the reference-only control experiment.
http:///www..oeiiiini.r.ffffieir stable oa
tty
33

---

Published as a conference paper at ICLR 2024
Number steps
110
steps PLMS4
(a) PLMS4
(b) with IIB 3.8
(c) GHVB3.1
Prompt: "a close up image of (s1) in conference"
Number steps
30
110
150
190
1,000 steps PLMS4
(a) PLMS4
(b) with HB 3.8
(c) GHVV3.1
(b) Prompt: "a close up image of (s1) under a tree"
Number of steps
30
70
110
150
190
1,000 steps PLMS4
(a) PLMSS4
(b)
(O) GHVB3.1
(c) Prompt: "a close up image (sl) in kitchen"
Figure 23: Comparison of samples generated from fine-tuned Stable Diffusion model called DreamShaper with (a) PLMS4 (Liu et al., 2022a). (b) PLMS4 with HB 0.8, and (c) GHVB3.1 using different number of sampling steps and guidance scale 15.
34

---

Published as a conference paper at ICLR 2024
Number steps
20
30
110
1,000 steps PLMS4
(a) -
(b) with HB 0.9
(a) Prompt: "a close up image (s1) in conference"
Number steps
110
1,000 steps
(a) HPM--Slvert+
(b) with HB 0.9
(b) Prompt: "a close up image (s1) under a tree"
Number steps
50
110
1,000 steps PLMS4
(a) DPM Solvert+
(b) with HB
(c) Prompt: "a close image of (s1) in forest"
Figure 24: Comparison of samples generated from fine-tuned Stable Diffusion model called DreamShaper with (a) DPM-Solver++ (Lu et al. 2022b) and (b) DPM-Solver++ with HB 0.9 using different number of sampling steps and guidance scale $ 15.
35

---

Published as a conference paper at ICLR 2024
Number of steps
30
70
110
150
190
Reference image
(a) PLMS4
1,000 steps PT MS4
(b) with HB 0.8
Number staps
70
110
Reference image
( PLMS4
steps PLMS4
(b)
Number steps
70
110
150
190
Reference image
(a) PLMS4
1,000 steps PLMS4
(c) GHVB3.1
Number of steps
30
70
110
150
190
Reference image
PLMS4
1,000 steps PLAMA4
GHVBB.1
Figure 25: Comparison samples generated from fine-tuned Stable Diffusion model named Any-
thing V3.0 in reference-only control pipeline with (a) PLMS4 (Liu et al.) 2022a), (b) PLMS4 with HB 0.8, and (c) GHVB3.1 using different number of sampling steps and guidance scale 12. Reference images are available on Hugging Face [8]
Prompt: "Igir!"
36

---

Published as a conference paper at ICLR 2024
P FACTORS CONTRIBUTING TO ARTIFACT OCCURRENCE IN FINE-TUNED DIFFUSION MODELS
We investigate three factors that influence the occurrence of divergence artifacts in diffusion sam- pling: number of steps, guidance scale, and choice diffusion models. We present results obtained from Stable Diffusion 1.5 (Figure and other models fine-tuned generate images in specific styles: Midjourney (Figure 29 Japanese animation (Figure 311 and photorealistic (Figure
Our observations reveal that insufficient numbers of steps and high guidance scales positively cor- relate with the presence of divergence artifacts in the generated samples. Fine-tuned models exhibit higher sensitivity these factors, resulting in greater incidence of artifacts compared to Stable Diffusion 1.5. Consistent with the findings presented in Section [.2] reducing the number of steps increases the likelihood of artifact occurrence. Furthermore, increasing the guidance scale amplifies the magnitude of eigenvalues, contributing to the presence of artifacts.
The choice of diffusion model also impacts artifact occurrence. Stable Diffusion 1.5 produces the fewest artifacts compared to the fine-tuned models. Among these fine-tuned models, Openjourney, which yields results similar to those obtained from Stable Diffusion 1.5, exhibits the lowest occur- rence of artifacts. This suggests that extensive changes to the model may alter the eigenvalues and lead to an increased presence of artifacts.
Additionally, present the results of our techniques for mitigating divergence artifacts in Figure 28 30 32 and 34 The parameter that offers the optimal tradcoff between reducing artifacts and mantaining accuracy depends on the chosen guidance scale and diffusion model.
Q ABLATION STUDY ON HB MOMENTUM
This section provides a comprehensive analysis of the convergence speed of Polyak's Heavy Ball (HB) momentum, which can be incorporated into existing diffusion sampling methods by modifying a few lines of code. Specifically, we compare the results obtained from several methods with and without HB momentum with the targets generated from the I, 1,000-step DDIM, using LPIPS in the image space and L2 in the latent space. We then estimate their orders of convergence, as explained in Appendix]] ! The results are visually presented in Figure 26
In contrast to the interpolation-lik behavior observed in Figure 12 for GHVB, observe that the use of HB momentum leads to an increase both the LPIPS score and the L2 distance when selecting values of B that less than 1. This is when 33 below 0.7. These findings indicate a deviation from the desired convergence behavior, highlighting potential decrease in solution accuracy, even though HB momentum has been shown to successfully mitigate divergence artifacts.
than the 1dtorder method DDIM
Additionally, we find that the numerical orders of convergence also tend approach the value. These observations align with our analysis in Theorem [ of Appendix H indicating that when B deviates from 1, the employed approach exhibits 1*-order convergence and is unable achieve high-order convergence. These conclusions emphasize the importance of carefully considering the choice of in order to strike balance between convergence speed and solution quality. Further details and insights into the performance of the HB momentum approach can be obtained from Figure 26] enhancing our understanding of behavior within the context of studied problem.
And
1190 M5.3
Figure 26: Comparison of LPIPS, mean L.2 distance, and order of convergence of HB when using different damping coefficients. Statistical means are averaged from 160 initial latent codes.
37

---

Published as a conference paper at ICLR 2024
"A post-apocalyptic world with ruined buildings, overgrown vegetation, and red sky"
"A girl standing in park in Japanese animation style"
22.5
22.5
10 steps
steps
20 steps
steps
Figure 27: Comparison of samples generated from Stable Diffusion 1.5 2022a) with different sampling steps and guidance scale.
19
Juing PLMS4 (Liu al,
apocaapoci world ruined buildings, overgrown vegetation, and sky"
girl standing park Jupanese unimation style"
7.5
15
22.5
7.5
15
22.5
10
15 steps
steps
40 steps
Figure 28: Comparison samples generated from Stable Diffusion 1.5 using PLMS4 with HB under various sampling steps and guidance scale Specifically, we employ 0.9 for = 7.5, 3 = 0.8 for s = 15, and 3 0.6 for & = 22.5 to account for the varying degrees of artifact manifestation associated with each guidance scale.
"ttps:://ogginniooo...///mmm//// -ww.iiii----fuooo------oo Ffu
38

---

Published as a conference paper at ICLR 2024
"A post-apocalyytic world with ruined buildings, overgrown vegetation, and red sky"
"A girl standing in park in Japanese animation style"
22.5
22.5
10 steps
steps
20 steps
steps
Figure 29: Comparison of samples generated from Openjourney usiin PLMS4 (Liu al. 2022a with different sampling steps and guidance scale.
apocalpoii world ruined buildings, overgrown vegetation, and sky"
girl standing park Jupanese unimation style"
7.5
= 15
22.5
7.5
15
22.5
10
15 steps
steps
40 steps
Figure 30: Comparison of samples generated from Openjourney using PLMS4 with HB 8 under various sampling steps and guidance scale Specifically, employ 3 0.8 for 7.5, B = 0.6 for s = 15, and / 0.6 for 22.5 to account for the varying degrees artifact manifestation associated with each guidance scale.
39

---

Published as a conference paper at ICLR 2024
"A post-ppocolyptic world with ruined buildings, overgrown vegetation, and red sky"
"A girl standing in park in Japanese animation style"
22.5
22.5
10 steps
steps
20 steps
steps
Figure 31: Comparison samples generated from Waifu Diffusion V1.4 IT using PLMS4 (Liu al. 2022a) with different sampling steps and guidance scale.
"A post-apocalypiic world ruined buildings, overgrown vegetation, and sky"
girl standing park Jupanese animation style"
7.5
15
22.5
7.5
15
22.5
10
15 steps
steps
40 steps
Figure 32: Comparison of samples generated from Waifu Diffusion VL.4 using PLMS4 with HB under various sampling steps and guidance scale Specifically, we employ = 0.8 for = 7.5, = 0.7 for s 15, and 3 0.6 for & = 22.5 to account for the varying degrees of artifact manifestation associated with each guidance scale.
httpp:///oggiingaee..c//aaaaeae fu-di ffusion
40

---

Published as a conference paper at ICLR 2024
"A post-apocalyptic world with ruined buildings, overgrown vegetation, and red sky"
"A girl standing in park in
Japanese animation style"
7.5
22.5
10 steps
steps
20 steps
steps
Figure 33: Comparison of samples generated from Dreamlike Photoreal 2.0 ct al. 2022a) with different sampling steps and guidance scale.
uusing PLMS4 (Liu
"A post upocalyptic world with ruined buildings, overgrown vegetation, and red sky"
"A girl standing in park in Japanese animation style"
10 steps
20 steps
steps
Figure 34: Comparison of samples generated from Dreamlike Photoreal V2.0 using PLMS4 with HB 3 under various sampling steps and guidance scale s. Specifically, we employ 3 0.7 for 7.5, 3 0.6 for 15, and 0.6 for 22.5 to account for the varying degrees of artifact manifestation associated with each guidance scale.
httpss://uuggiinffcccccc////// ike-art/dream ke peepphotoeell--e -2.
41

---

Published as a conference paper at ICLR 2024
"A post-apocolyptic world with ruined buildings, overgrown vegetation, and red sky" 7.5
"A girl standing in park in
Japanese animation style"
15
22.5
22.5
10 steps
steps
20 steps
steps
Figure 35: Comparison of samples generated from Stable Diffusion 1.5 using GHVB(3.0+ 9) under
various sampling steps and guidance scale s. Specifically, we employ 3 = 0.5 for s = 7.5, 3 = 0.2 for s 15, and 3 0.1 for 22.5 to account for the varying degrees of artifact manifestation associated with each guidance scale.
R ABLATION STUDY ON NESTEROV MOMENTUM
In Appendix G we investigated the potential of incorporating different types of momentum, such as Nesterov's ' momentum, into existing diffusion sampling methods to mitigate divergence artifacts. Similar the analysis conducted Section 5.4 and Appendix Q the primary objective of this section is to explore the convergence speed of Nesterov's momentum by comparing two key metries: LPIPS in the image space and L2 in the latent space.
Figure presents the results, which reveal intriguing parallels with the behavior of HB momentum observed in Figure 26 When Nesterov's 's momentum is applied to the PLMS2 method, the solution accuracy progressively diminishes the value of f deviates from as indicated by the correspond- ing increase in both LPIPS and 1.2 metrics. Notably, the model's accuracy drops below that of the DDIM when falls below 0.5. Qualitative comparisons of these two types of momentum shown in Figure 37
Furthermore, our analysis of the order of convergence demonstrates that Nesterov's momentum does not achieve a high order of convergence, similar to HB momentum. These findings emphasize the importance of carefully considering the choice of momentum method, along with the specific values assigned to 3, in order to strike an optimal balance between convergence speed and solution quality.
ABLATION STUDY ON MAGNITUDE SCORE
In this section. objective is to provide further verification and justification of the experiment ducted in Section by exploring various parameter settings for the magnitude score and assessing their effects on the selected model.
42

---

Published as a conference paper at ICLR 2024
A380 ardl
RAD CTCC
Figure 36: Comparison of LPIPS, mean L2 distance, and order of convergence of Nesterov's mo- mentum when using different damping coefficients. Statistical means are averaged from 160 initial latent codes.
S.1 RESULTS WITH ALTERNATIVE PARAMETER SETTINGS
To gain deeper insights into the integration of momentum into sampling methods, we analyze the results of the magnitude scores depicted in Figure 7 (Section 5.1). This analysis involves varying the threshold r and the kernel size for max- poxoppooiig the calculation of the magnitude score. By investigating different parameter settings, we aim to validate the outcomes of the experiment and uncover the scaling impact of the magnitude score. The results, shown in Figure 38 highlight influence of threshold and kernel size k on the magnitude score. It is important to note that while extreme values of may introduce ambiguity interpreting the outcomes, the overall observed trends remain consistent.
S.2 RESULTS ON ALTERNATIVE MODELS
In this section, we present the findings from our analysis conducted on alternative diffusion models, namely Stable Diffusion 1.5, Waifu Diffusion V1.4, and Dreamlike Photoreal 2.0. The primary aim of this investigation is to assess the impact of different models on the magnitude score and determine whether the trends identified in Section 5.1 hold across diverse model architectures.
For this analysis, we employed the same magnitude score parameters as in Section 5.1 of our examination are illustrated in Figure 339 which showcases the magnitude scores for each model. One important observation is that the change in model architecture only affects the scale of magnitude score, while the overall trend remains consistent across all models.
The results
STATISTICAL REPORTS
In this section, we present detailed statistical reports for the experiments conducted Section 5.1 and Section 4.2 These reports provide detailed information, including mean values and their corre- sponding 95% confidence intervals, to offer a thorough understanding the experimental results.
Firstly, focus on the experiment related to mitigating the magnitude score in Section 5.1 The results depicted Figure are presented in Table[8] Additionally, the outcomes illustrated in Figure
8 are reported in Table 9 For the ablation study of GHVB in Section 4.2 provide the results
shown in Figure [ in Table I0 and report the findings in Figure 12 in Table
Furthermore, include a runtime comparison of each sampling method in Table 12 detailing the wall clock time required for each method. The results indicate that all the methods exhibit similar sampling times, ensuring fair comparison across the different approaches.
43

---

Published as a conference paper at ICLR 2024
8 0.2
B = 0.4
B 0.6
B == 0.8
8 1.0
HB
8 0.2
B 0.4
3= 0.6
9 == 0.8
(b) w/ NT
(a) Prompt: beautiful illustration a schoolgirl riding her bicycle to school in small village"
0.2
0.4
0.6
0.8
3 1.0
0.2
0.4
8 0.6
# 0.8
NT
(b) Prompt:
"A beautiful illustration of a cozy cafe in a futuristic city"
8 0.2
0.4
3 U.6
3 0.8
9 1.0
(a) w/ HB
3 0.2
0 0.4
3
9 0.8
(b) NT
(c) Prompt: "A painting of an old European castle in deep forest with Blood Moon in the background"
Figure 37: Comparison between two variations momentum: (a) Polyak's Heavy Ball (HB) and (b) Nesterov (NT). These momentum variations are applied to PLMS4 (Liu et al. 2022a)
fine- tuned Stable Diffusion model called Anything V4 with 15 sampling steps and a guidance scale of 15. Both variations effectively reduce artifacts. However, the choice of the effectiveness parameter
3 might differ due to the distinct shapes of their respective stability regions.
44

---

Published as a conference paper at ICLR 2024
Threshold
1.5
3.0
10.0
Fiets
Fiets
kernel size
f:
kernel
1: ===4 4
Fuld
kermel size
1:
Figure 38: Comparison of magnitude scores on Anything V4 on different combinations of threshold and kernel size k used in max-pooling. #samples 160
Stable Diffusion 1.5
Waifu Diffusion V1.4
Dreamlike Photoreal 2.0
2.4"
PUM64 WIF 0.3"
FUMSA CD!
FUNS4 WiT 2.S' FUND4 win c. R.NSA
2.2
PUM64
FLPS4 H3
2.50-
2225
belva
2.8
210
125-
0.5
0.5-
20
10
Namber
Mamrar Stigo
NAmmoer Stees
PLMSA
PJM44
2.00
GW495.7
1.50
GH68231
3.50-
GAM83.7
May May-fudessos
GHV52.7
1.2
DOM
1:25
GH1/5.77
BG/M
1.00
1.00
0.79
2.0
0.50
20
22
80
35
40
Murbbrr Steps
Numter of Stude
Numeer of Grope
Figure 39: Comparison of magnitude scores in different diffusion models. #samples = 160
Number of steps
Method
10
255
30
DPM
1.113 + .090
1.369
+ .102
1.087 + .083 0.972 + .075 0.919 + .068 0.869 + .067
DPM w/ HB 0.8
0.974 + .078
1.057
.079 0.916
0.857 + -070 0.831 + .067 0.834 + .065
DPM w/ HB 0.9
1.043 + .082
1.186
+ .088 0.986 + .075 0.921
.073 0.867 + .068 0.844 + .065
PLMS4 w/ HB 0.8
1.958 .16 1.469
1 .105
1.213
.097
1.060
+ .087 0.963 + .076 0.838 1 .063
PLMS4 w/ HB 0.9
2.499
= 1.888
- .112
1.534
.112
1.270
1.116
.091 0.887
1 .066
PLMS4
3.149 + .16 2.460
.116 1.911
1 .115 1.597 .11 1.372 1 .106 0.957 1 .075
Table 8: 95% confidence intervals for the magnitude scores of HB (Figure) 7
45

---

Published as a conference paper at ICLR 2024
Number of steps
Method
10
15
20
25
30
60
DDIM
0.844 - .076 0.765 +
064
0.728
.60 0.744 + 063
0.761 + 062
0.778 + .62
GHVB2.1
1,238 = .097
0,924
0072
0.832 + .067 0,820 + 066
0.825 = 066
0,829 - .644
GHVB2.3
1,291 1991:.10 .10 0.952
0.872 + .070
0.836 - 064 0.831 = .66
0.828 = .62
GHVB2.5
1.392 + .103 1.016
0.907
.071 0.851 .669 0.842 + .065
0.845 .063
GHVB2.7
1.514 + .105
1.095
0.968
.077
0.877 + .68
0.864 + 067
+ 0.826
GHVB2.9
1.673 + .107
1.203 + 090
1.023 + .079
0.934 + .075
0.901
+ .071
0.835 + .66
PLMS4
3.149 + .116 2,460 1 .166 1991 + .115 1.597 + .116 1.372
+ .166 0.957 + .075
Table 9: 95% confidence intervals for the magnitude scores of GHVB (Figure
Number of steps
Method
10
20
40
80
160
320
640
DDIM
0.584 + .034 0.409 + .29
0.304 + .029 0.210 + .026 0.139 + .022 0.085 + 001 0.048 + .009
GHVB1.1
0.592
.034 0.406 + .29
0.295 + .00 0.189
.026 0.113
.220
0.054 + .010
0.019
.005
GHVB1.3
0.609
.035 0.410 .099
0.276 .22 0.158
023 0.086 + .017
0.030
- .007 0.009
.003
GHVB1.5
0.624 + .036 0.409 1 .029 0.261 1: .02 0.145 1 .20
0.067 1 .014 0.021
+ .000
0.006
1 .002
G.VV1.77
0.645 1: .037 0.411 .03
0.254
:228 0.133
1: .02
0.053
.011 0,016
1 .00
0,004 1 .002
GHVB1.9
0.663
.037 0.414
.030
0.246
0.123
.021
0.044
1 .009 0.013
1 .004
0.003
.001
PLMS2
0.676 .038 0.418 .030 0.246 028 0.119
021
0.041
.009
0.011
0.003
.001
Table 10: 95% confidence intervals for L2 norm of GHVB (Figrue |1
Number steps (kmew)
Method
40
80
160
320
640
GHVBO.5
0,247
.030
0.235
0.351 1: .045
0.450
.057
0.474
.057
GHVB1.5
0,550
_072
0.717 + .086
0,922 + .089
1.337 + .002
1.519 1002
GHVB2.5 0.624 + .077
1.121 .155
1.546
+ .132
1.906
+ .153
1.846 + .117
GHVB3.5 0.459 + .063 0.920 .00 1.877 = .10 1.960
+ .147
1.779 166
Table 11: 95% confidence intervals for the numerical orders of convergence of GHVB (Figure 12
Number of steps
Method
15
30
60
DPM-Solver++
2.49
4.84
9.54
DPM-Solver++ w/ HB 0.9
2.49
4.84
9.54
PLMS4
2.49
4.84
9.54
PLMS4 w/ HB 0.9
2.46
4.79
9.43
PLMS4 w/ NT 0.9
2.53
4.93
9.70
GHVB3.9
2.50
4.84
9.54
Table 12: Comparison of the average sampling time per image (in seconds) when using different
numbers of steps in Stable Diffusion 1.5 on single NVIDIA GeForce RTX 3080. The time differ- ences are marginal.
46


Reviews:
[{"review_id": "Ki346JFC86", "reviewer": "Reviewer_jPhx", "summary": "This submission suggests to use higher order numerical scheme (heavy ball momentum coupled with higher order multi-step methods in numerical ODE) to compute the diffusion process in computer vision.", "questions": "- Are metrics \"high-frequency error norm (HFEN)\" [MR image reconstruction from highly undersampled k-space data by dictionary learning, Ravishankar and Bresler, 2011] and Structural Similarity Index (SSIM) potentially relevant to measure the divergence artifacts (section 5.1)? If yes, then reporting evaluation results in these two metrics can be helpful.\n\n- Conceptually, I would like to understand better what authors mean by \"classifier-guided diffusion sampling\". What is the difference (conceptually and when it comes to implementation) between classifier-guidance and text-prompt based generation?", "limitations": "", "rating": 8, "confidence": 2, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "Authors' effort in experiments seem to be solid and thorough. \nAuthors have also been patient to review basics of stability concept in numerical ODEs.", "weaknesses": "I recommend that authors add a paragraph explaining what \"sampling\" means in the context of diffusion in the appendix, so that the content can be more self-contained. From what I understand about the main text, authors mean generating/inferring an image with trained diffusion models. This is not equivalent to the meaning of illustrating the distribution of all potentially generated images given underlying diffusion models.\n\n\nI also suggest that authors make a table to list all used numerical formats, explicitly, either in the main text or in appendix, to generate images. In this way, readers can associate the listed methods in each table/figure with specific algorithms. \nThe current presentation stops at a conceptual derivation of discrete update format instead of concrete update formula. In a similar spirit, it will be also helpful for authors to detail the setup of the training paradigm (specifically, what the loss function is for training).", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "r1ZYJrb5kX", "reviewer": "Reviewer_DpnC", "summary": "This paper considers the artifacts problem of ODE/SDE-based diffusion sampling. Authors thought that the divergence artifacts are caused by the stability regions of high-order numerical methods for solving ODEs and proposed two solutions for expanding the stability regions of current diffusion numerical methods, called Heavy Ball (HB) momentum and Generalized Heavy Ball. And in the case of low-step sampling, the proposed methods are effective in reducing artifacts. But the actual improvement on diffusion sampling acceleration is unlear.\n\n-------------\nPost-rebuttal: I read the rebuttal and thanks for the authors' efforts. I would like to keep my score.", "questions": "1.\tThe main difference between HB and GHVB is that HB calculates the moving average after summing high-order coefficients, whereas GHVB calculates it before the summation. Why does such a difference lead a larger stability region? \n2.\tCan the divergence artifacts be solved or mitigated by improving the dynamic range of pixel?\n3.\tWith additional training, what is the proposed methods’ complexity or cost?", "limitations": "", "rating": 6, "confidence": 2, "soundness": 3, "presentation": 2, "contribution": 3, "strengths": "1.\tThe divergence artifacts problem is theoretically linked with the stability region of high-order numerical solvers for ODEs. The insight is very helpful for the design of diffusion sampling methods. \n2.\tTo enlarge the stability region, authors proposed Heavy Ball (HB) and generalized Heavy Ball (GHVB) as two solution without any training. Experiments show that the divergence artifacts are great mitigated in a low-step sampling case.\n3.\tThis paper is well organized and solid in theory.", "weaknesses": "1.\tThe proposed method should be compared with the state-of-art methods in reducing divergence artifacts if it is a big challenge in diffusion models.\n2.\tThe stated motivation is diffusion model acceleration. Experiments are limited in comparing the results of few-step sampling, lacking clear numerical experiments in model acceleration.\n3.\tIt seems that the proposed methods show superior performance only in extremely low sampling steps. In the case of decent image quality, the improvement on sampling step is unclear.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "tYrh1K9WGE", "reviewer": "Reviewer_Z3BT", "summary": "This paper focuses on accelerating general diffusion sampling, where both unconditional and guided sampling are considered. Motivated by the observation that recent higher-order numerical methods would lead to diverging artifacts at lower sampling steps, the authors propose to incorporate heavy ball (HB) momentum into existing diffusion ODE solvers such as DPM++ and PLMS to mitigating their artifacts. In addition, an improved high-order version, namely generalized heavy Ball (GHVB) is also presented in this paper.  Experimental results have shown the effectiveness of this proposal.", "questions": "1), While the authors mentioned that the problem setup is more challenging in this paper than previous works, it is unclear what the challenges are. More discussions about why PLMS and DPM-Solver ++ perform worth than their original claims would strengthen this proposal.\n\n2), Given that the 1000-Step DDIM’s sample is considered the benchmark, it would be reasonable to include evaluation metrics such as L2, LIPIS, and FID comparing HB/GHVB to DDIM, as depicted in Figure 11.\n\n3), In Figure. 12, the authors attribute the inconsistency of GHVB 2.5 and 3.5 to estimated error or other sources of error without further justifications. It would be helpful to discuss this more for better understanding. \n\n4), Seems the comparisons and discussions between HB and GHVB are not sufficient in the paper’s current state. There is no clear cut which method is better for both conditional and unconditional diffusion sampling.", "limitations": "", "rating": 5, "confidence": 3, "soundness": 3, "presentation": 2, "contribution": 3, "strengths": "1), Both pixel-based and latent diffusion models are considered in this paper.\n\n2), The presentation is overall easy to follow.\n\n3), Good practical extension to DPM++ and PLMS.\n\n4), The literature over existing high order ODE solvers seems up to date.", "weaknesses": "1), The technical novelty behind this work seems to be not significant. The main techniques used in this paper are directly borrowed from Polyak’s heavy ball (HB) momentum method, a conventional optimization algorithm. Besides, the main improvements of this work are built based on DPM++ and PLMS.\n\n2), While two methods are proposed in the same paper, it is unclear which one should be used under what circumstances. The paper only gives some vague statements without comprehensive comparison. \n\n3), While guided diffusion sampling is considered, the effectiveness of the HB/GHVB under different scaling factor “s” is not well discussed.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "iomCTP2IVu", "reviewer": "Reviewer_uXih", "summary": "The paper introduces the Heavy Ball (HB) momentum into diffusion numerical methods to expand the stability regions. Meanwhile, the authors propose the high-order method, Generalized Heavy Ball (GHVB), to select the suitable method. Experiments show that the proposed HB and GHVB improves existing on both pixel-based and latent-based diffusion model in reducing artifacts and improving image quality.", "questions": "1. How effective are the proposed methods in extremely small generation steps, such as one or two?", "limitations": "", "rating": 8, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "1. The authors introduce the Heavy Ball (HB) momentum into existing diffusion methods to expand the stability regions. And they propose a high-order method, Generalized Heavy Ball (GHVB), to trade off between accuracy and artifact suppression.\n2. The analyses are adequate. Through visualization and theoretical analysis, it is discovered that the small stability regions lead to model artifacts.\n3. The experiments are extensive. The authors apply HB and GHVB on pixel-based and latent-based diffusion models to prove the effectiveness of the proposed method.\n4. The authors also provide the code, which shows the solidness of the work.", "weaknesses": "The paper primarily experiments with 10 or more generation steps. But, it lacks analyses of extreme cases, such as one or two steps. It is suggested to evaluate the effectiveness of the proposed methods in these scenarios, e.g., one or two steps. For instance, the consistency model [1] performs well in one- and few-step generation. How effective is the method proposed in this paper compared with CM?\n\n[1] Consistency Models, Song, Yang and Dhariwal, Prafulla and Chen, Mark and Sutskever, Ilya", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}]

Author rebuttal:

Decision: Accept (poster)