DPM-Solver: A Fast ODE Solver for Diffusion Probabilistic Model Sampling in Around 10 Steps

Introduction

Diffusion probabilistic models (DPMs) are emerging powerful generative models with promising performance on many tasks, such as image generation , video generation , text-to-image generation , speech synthesis and lossless compression . DPMs are defined by discrete-time random processes or continuous-time stochastic differential equations (SDEs) , which learn to gradually remove the noise added to the data points. Compared with the widely-used generative adversarial networks (GANs) and variational auto-encoders (VAEs) , DPMs can not only compute exact likelihood , but also achieve even better sample quality for image generation . However, to obtain high-quality samples, DPMs usually need hundreds or thousands of sequential steps of large neural network evaluations, thereby resulting in a much slower sampling speed than the single-step GANs or VAEs. Such inefficiency is becoming a critical bottleneck for the adoption of DPMs in downstream tasks, leading to an urgent request to design fast samplers for DPMs.

Existing fast samplers for DPMs can be divided into two categories. The first category includes knowledge distillation and noise level or sample trajectory learning . Such methods require a possibly expensive training stage before they can be used for efficient sampling. Furthermore, their applicability and flexibility might be limited. It might require nontrivial effort to adapt the method to different models, datasets, and number of sampling steps. The second category consists of training-free samplers, which are suitable for all pre-trained DPMs in a simple plug-and-play manner. Training-free samplers include adopting implicit or analytical generation process, advanced differential equation (DE) solvers and dynamic programming . However, these methods still require $\sim$ 50 function evaluations to generate high-quality samples (comparable to those generated by plain samplers in about 1000 function evaluations), thereby are still time-consuming.

In this work, we bring the efficiency of training-free samplers to a new level to produce high-quality samples in the “few-step sampling” regime, where the sampling can be done within around 10 steps of sequential function evaluations. We tackle the alternative problem of sampling from DPMs as solving the corresponding diffusion ordinary differential equations (ODEs) of DPMs, and carefully examine the structure of diffusion ODEs. Diffusion ODEs have a semi-linear structure — they consist of a linear function of the data variable and a nonlinear function parameterized by neural networks. Such structure is omitted in previous training-free samplers , which directly use black-box DE solvers. To utilize the semi-linear structure, we derive an exact formulation of the solutions of diffusion ODEs by analytically computing the linear part of the solutions, avoiding the corresponding discretization error. Furthermore, by applying change-of-variable, the solutions can be equivalently simplified to an exponentially weighted integral of the neural network. Such integral is very special and can be efficiently approximated by the numerical methods for exponential integrators .

Based on our formulation of solutions, we propose DPM-Solver, a fast dedicated solver for diffusion ODEs by approximating the above integral. Specifically, we propose first-order, second-order and third-order versions of DPM-Solver with convergence order guarantees. We further propose an adaptive step size schedule for DPM-Solver. In general, DPM-Solver is applicable to both continuous-time and discrete-time DPMs, and also conditional sampling with classifier guidance . Fig. 1 demonstrates the speedup performance of a Denoising Diffusion Implicit Models (DDIM) baseline and DPM-Solver, which shows that DPM-Solver can generate high-quality samples with as few as 10 function evaluations and is much faster than DDIM on the ImageNet 256x256 dataset . Our additional experimental results show that DPM-Solver can greatly improve the sampling speed of both discrete-time and continuous-time DPMs, and it can achieve excellent sample quality in around 10 function evaluations, which is much faster than all previous training-free samplers of DPMs.

Diffusion Probabilistic Models

We review diffusion probabilistic models and their associated differential equations in this section.

Under some regularity conditions, Song et al. show that the forward process in Eq. (2.2) has an equivalent reverse process from time $T$ to , starting with the marginal distribution $q_{T}(\bm{x}_{T})$:

where $\bar{\bm{w}}_{t}$ is a standard Wiener process in the reverse time. The only unknown term in Eq. (2.4) is the score function $\nabla_{\bm{x}}\log q_{t}(\bm{x}_{t})$ at each time $t$. In practice, DPMs use a neural network $\bm{\epsilon}_{\theta}(\bm{x}_{t},t)$ parameterized by $\theta$ to estimate the scaled score function: $-\sigma_{t}\nabla_{\bm{x}}\log q_{t}(\bm{x}_{t})$. The parameter $\theta$ is optimized by minimizing the following objective :

Samples can be generated from DPMs by solving the diffusion SDE in Eq. (2.5) with numerical solvers, which discretize the SDE from $T$ to . Song et al. proved that the traditional ancestral sampling method for DPMs can be viewed as a first-order SDE solver for Eq. (2.5). However, these first-order methods usually need hundreds of or thousands of function evaluations to converge , leading to extremely slow sampling speed.

2 Diffusion (Probability Flow) ODEs

When discretizing SDEs, the step size is limited by the randomness of the Wiener process [27, Chap. 11]. A large step size (small number of steps) often causes non-convergence, especially in high dimensional spaces. For faster sampling, one can consider the associated probability flow ODE , which has the same marginal distribution at each time $t$ as that of the SDE. Specifically, for DPMs, Song et al. proved that the probability flow ODE of Eq. (2.4) is

where the marginal distribution of $\bm{x}_{t}$ is also $q_{t}(\bm{x}_{t})$. By replacing the score function with the noise prediction model, Song et al. defined the following parameterized ODE (diffusion ODE):

Samples can be drawn by solving the ODE from $T$ to . Comparing with SDEs, ODEs can be solved with larger step sizes as they have no randomness. Furthermore, we can take advantage of efficient numerical ODE solvers to accelerate the sampling. Song et al. used the RK45 ODE solver for the diffusion ODEs, which generates samples in $\sim$ 60 function evaluations to reach comparable quality with a 1000-step SDE solver for Eq. (2.5) on the CIFAR-10 dataset . However, existing general-purpose ODE solvers still cannot generate satisfactory samples in the few-step ($\sim$ 10 steps) sampling regime. To the best of our knowledge, there is still a lack of training-free samplers for DPMs in the few-step sampling regime, and the sampling speed of DPMs is still a critical issue.

Customized Fast Solvers for Diffusion ODEs

As highlighted in Sec. 2.2, discretizing SDEs is generally difficult in high dimensions [27, Chap. 11] and it is hard to converge within few steps. In contrast, ODEs are easier to solve, yielding a potential for fast samplers. However, as mentioned in Sec. 2.2, the general black-box ODE solver used in previous work empirically fails to converge in few steps. This motivates us to design a dedicated solver for diffusion ODEs to enable fast and high-quality few-step sampling. We start with a detailed investigation of the specific structure of diffusion ODEs.

The key insight of this work is that given an initial value $\bm{x}_{s}$ at time $s>0$, the solution $\bm{x}_{t}$ at each time $t<s$ of diffusion ODEs in Eq. (2.7) can be simplified into a very special exact formulation which can be efficiently approximated.

Our first key observation is that a part of the solution $\bm{x}_{t}$ can be exactly computed by considering the particular structure of diffusion ODEs. The r.h.s. of diffusion ODEs in Eq. (2.7) consists of two parts: the part $f(t)\bm{x}_{t}$ is a linear function of $\bm{x}_{t}$, and the other part $\frac{g^{2}(t)}{2\sigma_{t}}\bm{\epsilon}_{\theta}(\bm{x}_{t},t)$ is generally a nonlinear function of $\bm{x}_{t}$ because of the neural network $\bm{\epsilon}_{\theta}(\bm{x}_{t},t)$. This type of ODE is referred to as semi-linear ODE. The black-box ODE solvers adopted by previous work are ignorant of this semi-linear structure as they take the whole $\bm{h}_{\theta}(\bm{x}_{t},t)$ in Eq. (2.7) as the input, which causes discretization errors of both the linear and nonlinear term. We note that for semi-linear ODEs, the solution at time $t$ can be exactly formulated by the “variation of constants” formula :

This formulation decouples the linear part and the nonlinear part. In contrast to black-box ODE solvers, the linear part is now exactly computed, which eliminates the approximation error of the linear term. However, the integral of the nonlinear part is still complicated because it couples the coefficients about the noise schedule (i.e., $f(\tau),g(\tau),\sigma_{\tau}$) and the complex neural network $\bm{\epsilon}_{\theta}$, which is still hard to approximate.

Our second key observation is that the integral of the nonlinear part can be greatly simplified by introducing a special variable. Let $\lambda_{t}\coloneqq\log(\alpha_{t}/\sigma_{t})$ (one half of the log-SNR), then $\lambda_{t}$ is a strictly decreasing function of $t$ (due to the definition of DPMs as discussed in Sec. 2.1). We can rewrite $g(t)$ in Eq. (2.3) as

Given an initial value $\bm{x}_{s}$ at time $s>0$, the solution $\bm{x}_{t}$ at time $t\in[0,s]$ of diffusion ODEs in Eq. (2.7) is:

We call the integral $\int e^{-\lambda}\hat{\bm{\epsilon}}_{\theta}(\hat{\bm{x}}_{\lambda},\lambda)\differential\lambda$ the exponentially weighted integral of $\hat{\bm{\epsilon}}_{\theta}$, which is very special and highly related to the exponential integrators in the literature of ODE solvers . To the best of our knowledge, such formulation has not been revealed in prior work of diffusion models.

Eq. (3.4) provides a new perspective for approximating the solutions of diffusion ODEs. Specifically, given $\bm{x}_{s}$ at time $s$, According to Eq. (3.4), approximating the solution at time $t$ is equivalent to directly approximating the exponentially weighted integral of $\hat{\bm{\epsilon}}_{\theta}$ from $\lambda_{s}$ to $\lambda_{t}$, which avoids the error of the linear terms and is well-studied in the literature of exponential integrators . Based on this insight, we propose fast solvers for diffusion ODEs, as detailed in the following sections.

2 High-Order Solvers for Diffusion ODEs

In this section, we propose high-order solvers for diffusion ODEs with convergence order guarantee by leveraging our proposed solution formulation Eq. (3.4). The proposed solvers and analysis are highly motivated by the methods of exponential integrators in the ODE literature.

Substituting the above Taylor expansion into Eq. (3.5) yields

where the integral $\int e^{-\lambda}\frac{(\lambda-\lambda_{t_{i-1}})^{n}}{n!}\differential\lambda$ can be analytically computed by repeatedly applying $n$ times of integration-by-parts (see Appendix B.2). Therefore, to approximate $\bm{x}_{t_{i-1}\to t_{i}}$, we only need to approximate the $n$-th order total derivatives $\hat{\bm{\epsilon}}_{\theta}^{(n)}(\hat{\bm{x}}_{\lambda},\lambda)$ for $n\leq k-1$, which is a well-studied problem in the ODE literature . By dropping the $\mathcal{O}(h_{i}^{k+1})$ error term and approximating the first $(k-1)$-th total derivatives with the “stiff order conditions” , we can derive $k$-th-order ODE solvers for diffusion ODEs. We name such solvers as DPM-Solver overall, and DPM-Solver-$k$ for a specific order $k$. Here we take $k=1$ for demonstration. In this case, Eq. (3.6) becomes

By dropping the high-order error term $\mathcal{O}(h_{i}^{2})$, we can obtain an approximation for $\bm{x}_{t_{i-1}\to t_{i}}$. As $k=1$ here, we call this solver DPM-Solver-1, and the detailed algorithm is as following.

For $k\geq 2$, approximating the first $k$ terms of the Taylor expansion needs additional intermediate points between $t$ and $s$ . The derivation is more technical so we defer it to Appendix B. Below we propose algorithms for $k=2,3$ and name them as DPM-Solver-2 and DPM-Solver-3, respectively.

Here, $t_{\lambda}(\cdot)$ is the inverse function of $\lambda(t)$, which has an analytical formulation for the practical noise schedule used in , as shown in Appendix D. The chosen intermediate points are ($s_{i}$, $\bm{u}_{i}$) for DPM-Solver-2 and $(s_{2i-1},\bm{u}_{2i-1})$ and $(s_{2i},\bm{u}_{2i})$ for DPM-Solver-3. As shown in the algorithm, DPM-Solver-$k$ requires $k$ function evaluations per step for $k=1,2,3$. Despite the more expensive steps, higher-order solvers ($k=2,3$) are usually more efficient since they require much fewer steps to converge, due to their higher convergence order. We show that DPM-Solver-$k$ is $k$-th-order solver, as stated in the following theorem. The proof is in Appendix B.

Finally, solvers with $k\geq 4$ need much more intermediate points as shown by previous work for exponential integrators. Therefore, we only consider $k$ from $1$ to $3$ in this work, while leaving the solvers with higher $k$ for future study.

3 Step Size Schedule

The proposed solvers in Sec. 3.2 need to specify the time steps $\{t_{i}\}_{i=0}^{M}$ in advance. We propose two choices of the time step schedule. One choice is handcrafted, which is to uniformly split the interval $[\lambda_{T}$, $\lambda_{0}$], i.e. $\lambda_{t_{i}}=\lambda_{T}+\frac{i}{M}(\lambda_{0}-\lambda_{T})$, $i=0,\dots,M$. Note that this is different from previous work which chooses uniform steps for $t_{i}$. Empirically, DPM-Solver with uniform time steps $\lambda_{t_{i}}$ can already generate quite good samples in few steps, where results are listed in Appendix E. As the other choice, we propose an adaptive step size algorithm, which dynamically adjusts the step size by combining different orders of DPM-Solver. The adaptive algorithm is inspired by and we defer its implementation details to Appendix C.

For few-step sampling, we need to use up all the number of function evaluations (NFE). When the NFE is not divisible by $3$, we firstly apply DPM-Solver-3 as much as possible, and then add a single step of DPM-Solver-1 or DPM-Solver-2 (dependent on the reminder of $K$ divided by $3$), as detailed in Appendix D. In the subsequent experiments, we use such combination of solvers with the uniform step size schedule for NFE $\leq 20$, and otherwise the adaptive step size schedule.

4 Sampling from Discrete-Time DPMs

Comparison with Existing Fast Sampling Methods

Here, we discuss the relationship and highlight the difference between DPM-Solver and existing ODE-based fast sampling methods for DPMs. We further briefly discuss the advantage of training-free samplers over those training-based ones.

Although motivated by entirely different perspectives, we show that the updates of DPM-Solver-1 and Denoising Diffusion Implicit Models (DDIM) are identical. By the definition of $\lambda$, we have $\frac{\sigma_{t_{i-1}}}{\alpha_{t_{i-1}}}=e^{-\lambda_{t_{i-1}}}$ and $\frac{\sigma_{t_{i}}}{\alpha_{t_{i}}}=e^{-\lambda_{t_{i}}}$. Plugging these and $h_{i}=\lambda_{t_{i}}-\lambda_{t_{i-1}}$ to Eq. (4.1) results in exactly a step of DPM-Solver-1 in Eq. (3.7). However, the semi-linear ODE formulation of DPM-Solver allows for principled generalization to higher-order solvers and convergence order analysis.

Recent work also show that DDIM is a first-order discretization of diffusion ODEs by differentiating both sides of Eq. (4.1). However, they cannot explain the difference between DDIM and the first-order Euler discretization of diffusion ODEs. In contrast, by showing that DDIM is a special case of DPM-Solver, we reveal that DDIM makes full use of the semi-linearity of diffusion ODEs, which explains its superiority over traditional Euler methods.

2 Comparison with Traditional Runge-Kutta Methods

One can obtain a high-order solver by directly applying traditional explicit Runge-Kutta (RK) methods to the diffusion ODE in Eq. (2.7). Specifically, RK methods write the solution of Eq. (2.7) in the following integral form:

and use some intermediate time steps between $[t,s]$ and combine the evaluations of $\bm{h}_{\theta}$ at these time steps to approximate the whole integral. The approximation error of explicit RK methods depends on $\bm{h}_{\theta}$, which consists of the error corresponding to both the linear term $f(\tau)\bm{x}_{\tau}$ and the nonlinear noise prediction model $\bm{\epsilon}_{\theta}$. However, the error of the linear term may increase exponentially because the exact solution of the linear term has an exponential coefficient (as shown in Eq. (3.1)). There are many empirical evidence showing that directly using explicit RK methods for semi-linear ODEs may suffer from unstable numerical issues for large step size. We also demonstrate the empirical difference of the proposed DPM-Solver and the traditional explicit RK methods in Sec. 5.1, which shows that DPM-Solver have smaller discretization errors than the RK methods with the same order.

3 Training-based Fast Sampling Methods for DPMs

Samplers that need extra training or optimization include knowledge distillation , learning the noise level or variance , and learning the noise schedule or sample trajectory . Although the progressive distillation method can obtain a fast sampler within 4 steps, it needs further training costs and loses part of the information in the original DPM (e.g., after distillation, the noise prediction model cannot predict the noise (score function) at every time step between $[0,T]$). In contrast, training-free samplers can keep all the information of the original model, and thereby can be directly extended to the conditional sampling by combining the original model and an external classifier (e.g. see Appendix D for the conditional sampling with classifier guidance).

Beyond directly designing fast samplers for DPMs, several works also propose novel types of DPMs which supports faster sampling. For instance, defining a low-dimensional latent variable for DPMs ; designing special diffusion processes with bounded score functions ; combining GANs with the reverse process of DPMs . The proposed DPM-Solver may also be suitable for accelerating the sampling of these DPMs, and we leave them for future work.

Experiments

In this section, we show that as a training-free sampler, DPM-Solver can greatly speedup the sampling of existing pre-trained DPMs, including both continuous-time and discrete-time ones, with both linear noise schedule and cosine noise schedule . We vary different number of function evaluations (NFE) which is the number of calls to the noise prediction model $\bm{\epsilon}_{\theta}(\bm{x}_{t},t)$, and compare the sample quality between DPM-Solver and other methods. For each experiment, We draw 50K samples and use the widely adopted FID score to evaluate the sample quality, where lower FID usually implies better sample quality.

Unless explicitly mentioned, we always use the solver combination with the uniform step size schedule in Sec. 3.3 if the NFE budget is less than 20, and otherwise the DPM-Solver-3 with the adaptive step size schedule in Sec. 3.3. We refer to Appendix D for other implementation details of DPM-Solver and Appendix E for detailed settings.

We firstly compare DPM-Solver with other continuous-time sampling methods for DPMs. The compared methods include the Euler-Maruyama discretization for diffusion SDEs , the adaptive step size solver for diffusion SDEs and the RK methods for diffusion ODEs in Eq. (2.7). We compare these methods for sampling from a pre-trained continuous-time “VP deep” model on the CIFAR-10 dataset with the linear noise schedule.

Fig. 2(a) shows the efficiency of compared solvers. We use uniform time steps with 50, 200, 1000 NFE for the diffusion SDE with Euler discretization, and vary the tolerance hyperparameter for the adaptive step size SDE solver and RK45 ODE solver to control the NFE. DPM-Solver can generate good sample quality within around 10 NFE, while other solvers have large discretization error even in 50 NFE, which shows that DPM-Solver can achieve $\sim$5 speedup of the previous best solver. In particular, we achieve 4.70 FID with 10 NFE, 3.75 FID with 12 NFE, 3.24 FID with 15 NFE, and 2.87 FID with 20 NFE, which is the fastest sampler on CIFAR-10.

As an ablation study, we also compare the second-order and third-order DPM-Solver and RK methods, as shown in Table 1. We compare RK methods for diffusion ODEs w.r.t. both time $t$ in Eq. (2.7) and half-log-SNR $\lambda$ by applying change-of-variable (see detailed formulations in Appendix E.1). The results show that given the same NFE, the sample quality of DPM-Solver is consistently better than RK methods with the same order. The superior efficiency of DPM-Solver is particularly evident in the few-step regime under 15 NFE, where RK methods have rather large discretization errors. This is mainly because DPM-Solver analytically computes the linear term, avoiding the corresponding discretization error. Besides, the higher order DPM-Solver-3 converges faster than DPM-Solver-2, which matches the order analysis in Theorem 3.2.

2 Comparison with Discrete-Time Sampling Methods

We use the method in Sec. 3.4 for using DPM-Solver in discrete-time DPMs, and then compare DPM-Solver with other discrete-time training-free samplers, including DDPM , DDIM , Analytic-DDPM , Analytic-DDIM , PNDM , FastDPM and Itô-Taylor . We also compare with GGDM , which uses the same pre-trained model but needs further training for the sampling trajectory. We compare the sample quality by varying NFE from 10 to 1000.

Specifically, we use the discrete-time model trained by $L_{\text{simple}}$ in on the CIFAR-10 dataset with linear noise schedule; the discrete-time model in on CelebA 64x64 with linear noise schedule; the discrete-time model trained by $L_{\text{hybrid}}$ in on ImageNet 64x64 with cosine noise schedule; the discrete-time model with classifier guidance in on ImageNet 128x128 with linear noise schedule; the discrete-time model in on LSUN bedroom 256x256 with linear noise schedule. For the models trained on ImageNet, we only use their “mean” model and omit the “variance” model. As shown in Fig. 2, on all datasets, DPM-Solver can obtain reasonable samples within 12 steps (FID 4.65 on CIFAR-10, FID 3.71 on CelebA 64x64 and FID 19.97 on ImageNet 64x64, FID 4.08 on ImageNet 128x128), which is $4\sim 16\times$ faster than the previous fastest training-free sampler. DPM-Solver even outperforms GGDM, which requires additional training.

Conclusions

We tackle the problem of fast and training-free sampling from DPMs. We propose DPM-Solver, a fast dedicated training-free solver of diffusion ODEs for fast sampling of DPMs in around 10 steps of function evaluations. DPM-Solver leverages the semi-linearity of diffusion ODEs and it directly approximates a simplified formulation of exact solutions of diffusion ODEs, which consists of an exponentially weighted integral of the noise prediction model. Inspired by numerical methods for exponential integrators, we propose first-order, second-order and third-order DPM-Solver to approximate the exponentially weighted integral of noise prediction models with theoretical convergence guarantee. We propose both handcrafted and adaptive step size schedule, and apply DPM-Solver for both continuous-time and discrete-time DPMs. Our experimental results show that DPM-Solver can generate high-quality samples in around 10 function evaluations on various datasets, and it can achieve $4\sim 16\times$ speedup compared with previous state-of-the-art training-free samplers.

Limitations and broader impact Despite the promising speedup performance, DPM-Solver is designed for fast sampling, which may be not suitable for accelerating the likelihood evaluations of DPMs. Besides, compared to the commonly-used GANs, diffusion models with DPM-Solver are still not fast enough for real-time applications. In addition, like other deep generative models, DPMs may be used to generate adverse fake contents, and the proposed solver may further amplify the potential undesirable influence of deep generative models for malicious applications.

Acknowledgements

This work was supported by National Key Research and Development Project of China (No. 2021ZD0110502); NSF of China Projects (Nos. 62061136001, 61620106010, 62076145, U19B2034, U1811461, U19A2081, 6197222, 62106120); Beijing NSF Project (No. JQ19016); Beijing Outstanding Young Scientist Program NO. BJJWZYJH012019100020098; a grant from Tsinghua Institute for Guo Qiang; the NVIDIA NVAIL Program with GPU/DGX Acceleration; the High Performance Computing Center, Tsinghua University; the Fundamental Research Funds for the Central Universities, and the Research Funds of Renmin University of China (22XNKJ13). J.Z is also supported by the XPlorer Prize.

References

Checklist

Do the main claims made in the abstract and introduction accurately reflect the paper’s contributions and scope? [Yes]

Did you describe the limitations of your work? [Yes] See section 6.

Did you discuss any potential negative societal impacts of your work? [Yes] See section 6.

Have you read the ethics review guidelines and ensured that your paper conforms to them? [Yes]

If you are including theoretical results…

Did you state the full set of assumptions of all theoretical results? [Yes] See Appendix B.

Did you include complete proofs of all theoretical results? [Yes] See Appendix B.

Did you include the code, data, and instructions needed to reproduce the main experimental results (either in the supplemental material or as a URL)? [Yes] Code is attached in the supplemental materials, with the appendix.

Did you specify all the training details (e.g., data splits, hyperparameters, how they were chosen)? [Yes] Our method is training-free. But we also report the hyperparameters for evaluations used in our proposed solver.

Did you report error bars (e.g., with respect to the random seed after running experiments multiple times)? [No] We observe that the standard deviation of the FID evaluations of DPM-Solver are rather small (mainly less than 0.01) because the FID is already averaged over 50K samples, following existing work . The small standard deviation does not change the conclusion.

Did you include the total amount of compute and the type of resources used (e.g., type of GPUs, internal cluster, or cloud provider)? [Yes] The GPU type and amount is detailed in Appendix E.

If you are using existing assets (e.g., code, data, models) or curating/releasing new assets…

If your work uses existing assets, did you cite the creators? [Yes]

Did you mention the license of the assets? [Yes] See Appendix E

Did you include any new assets either in the supplemental material or as a URL? [Yes] We include our code in the supplemental materials.

Did you discuss whether and how consent was obtained from people whose data you’re using/curating? [No] All of the datasets used in the experiments are publicly available.

Did you discuss whether the data you are using/curating contains personally identifiable information or offensive content? [Yes] We mentioned the human privacy issues of the ImageNet dataset in Appendix E.

If you used crowdsourcing or conducted research with human subjects…

Did you include the full text of instructions given to participants and screenshots, if applicable? [N/A]

Did you describe any potential participant risks, with links to Institutional Review Board (IRB) approvals, if applicable? [N/A]

Did you include the estimated hourly wage paid to participants and the total amount spent on participant compensation? [N/A]

Appendix A Sampling with Invariance to the Noise Schedule

In this section, we discuss more about the exact solution in Proposition 3.1 and give some insights about the formulation. Below we firstly restate the proposition w.r.t. $\lambda$ (i.e. the half-logSNR).

Given an initial value $\hat{\bm{x}}_{\lambda_{s}}$ at time $s$ with the corresponding half-logSNR $\lambda_{s}$, the solution $\hat{\bm{x}}_{\lambda_{t}}$ at time $t$ of diffusion ODEs in Eq. (2.7) with the corresponding half-logSNR $\lambda_{t}$ is:

In the following subsections, we will show that such formulation decouples the model $\bm{\epsilon}_{\theta}$ from the specific noise schedule, and thus is invariant to the noise schedule. Moreover, such change-of-variable for $\lambda$ in Proposition 3.1 is highly related to the maximum likelihood training of diffusion models . We show that both the maximum likelihood training and the sampling of diffusion models have invariance formulations that are independent of the noise schedule.

In this section, we show that Proposition 3.1 can decouples the exact solutions of the diffusion ODEs from the specific noise schedules (i.e. choice of the functions $\alpha_{t}=\alpha(t)$ and $\sigma_{t}=\sigma(t)$). Namely, given a starting point $\lambda_{s}$, a ending point $\lambda_{t}$, an initial value $\hat{\bm{x}}_{\lambda_{s}}$ at $\lambda_{s}$ and a noise prediction model $\hat{\bm{\epsilon}}_{\theta}$, the solution of $\hat{\bm{x}}_{\lambda_{t}}$ is invariant of the noise schedule between $\lambda_{s}$ and $\lambda_{t}$.

We firstly consider the VP type diffusion models, which is equivalent to the original DDPM . For VP type diffusion models, we always have $\alpha_{t}^{2}+\sigma_{t}^{2}=1$, so defining the noise schedule is equivalent to defining the function $\alpha_{t}=\alpha(t)$ (For example, DDPM uses a noise schedule such that $\beta(t)=\frac{\differential\log\alpha_{t}}{\differential t}$ is a linear function of $t$, and i-DDPM uses a noise schedule such that $\beta(t)=\frac{\differential\log\alpha_{t}}{\differential t}$ is a cosine function of $t$). As $\lambda_{t}=\log\alpha_{t}-\log\sigma_{t}$, we have $\alpha_{t}=\sqrt{\frac{1}{1+e^{-2\lambda_{t}}}}$ and $\sigma_{t}=\sqrt{\frac{1}{1+e^{2\lambda_{t}}}}$. Thus, we can directly compute the $\alpha_{t}$ and $\sigma_{t}$ for a given $\lambda_{t}$. Denote $\hat{\alpha}_{\lambda}\coloneqq\sqrt{\frac{1}{1+e^{-2\lambda}}}$, we have

We should notice that the integrand $e^{-\lambda}\hat{\bm{\epsilon}}_{\theta}(\hat{\bm{x}}_{\lambda},\lambda)$ is a function of $\lambda$, so its integral from $\lambda_{s}$ to $\lambda_{t}$ is only dependent on the starting point $\lambda_{s}$, the ending point $\lambda_{t}$ and the function $\hat{\bm{\epsilon}}_{\theta}$, which is independent of the intermediate values. As other coefficients ($\hat{\alpha}_{\lambda_{s}}$ and $\hat{\alpha}_{\lambda_{t}}$) are also only dependent on the starting point $\lambda_{s}$ and the ending point $\lambda_{t}$, we can conclude that $\hat{\bm{x}}_{\lambda_{t}}$ is invariant of the specific choice of the noise schedules. Intuitively, this is because we converts the original integral of time $t$ in Eq. (3.1) to the integral of $\lambda$, and the functions $f(t)$ and $g(t)$ are converted to an analytical formulation $e^{-\lambda}$, which is invariant to the specific choices of $f(t)$ and $g(t)$. Finally, for other types of diffusion models (such as the VE type and the subVP type), they are all equivalent to the VP type by equivalently rescaling the noise prediction models, as proved in . Therefore, the solutions of these types also have such property.

In summary, Proposition 3.1 decouples the solution of diffusion ODEs from the noise schedules, which gives us an opportunity to design tailor-made samplers for DPMs. In fact, as shown in Sec. 3.2, the only approximation of the proposed DPM-Solver is about the Taylor expansion of the neural network $\hat{\bm{\epsilon}}_{\theta}$ w.r.t. $\lambda$, and DPM-Solver analytically computes other coefficients (which are corresponding to the specific noise schedules). Intuitively, DPM-Solver keeps the known information as much as possible, and only approximates the intractable integral of the neural network, so it can generate comparable samples within much fewer steps.

A.2 Choosing Time Steps for λ𝜆\lambda is Invariant to the Noise Schedule

As mentioned in Appendix A.1, the formulation of Proposition 3.1 decouples the sampling solution from the noise schedule. The solution depends on the starting point $\lambda_{s}$ and the ending point $\lambda_{t}$, and is invariant to the intermediate noise schedule. Similarly, the updating equations of the algorithm of DPM-Solver are also invariant to the intermediate noise schedule. Therefore, if we have chosen the time steps $\{\lambda_{i}\}_{i=0}^{M}$, then the solution of DPM-Solver is also determined and is invariant to the intermediate noise schedule.

A simple way for choosing time steps for $\lambda$ is uniformly splitting $[\lambda_{T},\lambda_{\epsilon}]$, which is the setting in our experiments. However, we believe that there exists more precise ways for choosing the time steps, and we leave it for future work.

A.3 Relationship with the Maximum Likelihood Training of Diffusion Models

Interestingly, the maximum likelihood training of diffusion SDEs in continuous time also has such invariance property . Below we briefly review the maximum likelihood training loss of diffusion SDEs, and then propose a new insight for understanding diffusion models.

Denote the data distribution as $q_{0}(\bm{x}_{0})$, the distribution of the forward process at each time $t$ as $q_{t}(\bm{x}_{t})$, the distribution of the reverse process at each time $t$ as $p_{t}(\bm{x}_{t})$ with $p_{T}=\mathcal{N}(\bm{0},\bm{I})$. In , it is proved that the KL-divergence between $q_{0}$ and $p_{0}$ can be bounded by a weighted score matching loss:

where $\bm{x}_{t}=\alpha_{t}\bm{x}_{0}+\sigma_{t}\bm{\epsilon}$ and $C$ is a constant independent of $\theta$. As shown in Sec. 3.1, we have

so by applying change-of-variable w.r.t. $\lambda$, we have

which is equivalent to the importance sampling trick in [41, Sec. 5.1] and the continuous-time diffusion loss in [10, Eq. (22)]. Compared to Proposition 3.1, we can find that the sampling and the maximum likelihood training of diffusion models can both be converted to an integral w.r.t. $\lambda$, such that the formulation is invariant to the specific noise schedules, and we summarize it in Table 2. Such invariance property for both training and sampling brings a new insight for understanding diffusion models. For instance, we can directly define the noise prediction model $\bm{\epsilon}_{\theta}$ w.r.t. the (half-)logSNR $\lambda$ instead of the time $t$, then the training and sampling for diffusion models can be done without further choosing any ad-hoc noise schedules. Such finding may unify the different ways of the training and the inference of diffusion models, and we leave it for future study.

Appendix B Proof of Theorem 3.2

Throughout this section, we denote $\bm{x}_{s}$ as the solution of the diffusion ODE Eq. (2.7) starting from $\bm{x}_{T}$. For DPM-Solver-$k$ we make the following assumptions:

The total derivatives $\frac{\differential^{j}\hat{\bm{\epsilon}}_{\theta}(\hat{\bm{x}}_{\lambda},\lambda)}{\differential\lambda^{j}}$ (as a function of $\lambda$) exist and are continuous for $0\leq j\leq k+1$.

The function $\bm{\epsilon}_{\theta}(\bm{x},s)$ is Lipschitz w.r.t. to its first parameter $\bm{x}$.

B.2 General Expansion of the Exponentially Weighted Integral

Firstly, we derive the Taylor expansion of the exponentially weighted integral. Let $t<s$ and then $\lambda_{t}>\lambda_{s}$. Denote $h\coloneqq\lambda_{t}-\lambda_{s}$, and the $k$-th order total derivative $\hat{\bm{\epsilon}}_{\theta}^{(k)}(\hat{\bm{x}}_{\lambda},\lambda)\coloneqq\frac{\differential^{k}\hat{\bm{\epsilon}}_{\theta}(\hat{\bm{x}}_{\lambda},\lambda)}{\differential\lambda^{k}}$. For $n\geq 0$, the $n$-th order Taylor expansion of $\hat{\bm{\epsilon}}_{\theta}(\hat{\bm{x}}_{\lambda},\lambda)$ w.r.t. $\lambda$ is

To expand the exponential integrator, we further define :

and it satisfies $\varphi_{k}(0)=\frac{1}{k!}$ and a recurrence relation $\varphi_{k+1}(z)=\frac{\varphi_{k}(z)-\varphi_{k}(0)}{z}$. By taking the Taylor expansion of $\hat{\bm{\epsilon}}_{\theta}(\hat{\bm{x}}_{\lambda},\lambda)$, the exponential integrator can be rewritten as

So the solution of $\bm{x}_{t}$ in Eq. (3.4) can be expanded as

Finally, we list the closed-forms of $\varphi_{k}$ for $k=1,2,3$:

B.3 Proof of Theorem 3.2 when k=1𝑘1k=1

Taking $n=0,t=t_{i},s=t_{i-1}$ in Eq. (B.4), we obtain

By Assumption B.2 and Eq. (3.7), it holds that

B.4 Proof of Theorem 3.2 when k=2𝑘2k=2

We prove the discretization error of the general form of DPM-Solver-2 in Algorithm 4.

First, we consider the following update for $0<t<s<T,h:=\lambda_{t}-\lambda_{s}$.

In this remaining part we prove that $\bar{\bm{x}}_{t}=\bm{x}_{t}+\mathcal{O}(h^{3})$.

Note that by the Lipschitzness of $\bm{\epsilon}_{\theta}$ w.r.t. $\bm{x}$ (Assumption B.2),

where the last equation follows from a similar argument in the proof of $k=1$. Since $e^{h}-1=\mathcal{O}(h)$, the second term of the above display is $\mathcal{O}(h^{3})$.

As $\lambda_{s_{1}}-\lambda_{s}=r_{1}h$, $\varphi_{i}(h)=(e^{h}-1)/h$ and $\varphi_{2}(h)=(e^{h}-h-1)/h^{2}$, we find

Then, the proof is completed by noticing that

B.5 Proof of Theorem 3.2 when k=3𝑘3k=3

As in Appendix B.4, it suffices to show that the following update has error $\bar{\bm{x}}_{t}=\bm{x}_{t}+\mathcal{O}(h^{4})$ for $0<t<s<T$ and $h=\lambda_{s}-\lambda_{t}$.

Similar to the proof in Appendix B.4, since $\frac{e^{r_{2}h-1}}{r_{2}h}-1=\mathcal{O}(h)$ and $\bar{\bm{u}}_{1}=\bm{x}_{s_{1}}+\mathcal{O}(h^{2})$, then

Let $h_{2}=r_{2}h$, then following the same line of arguments in the proof of Appendix B.4, it suffices to check that

which holds by applying Taylor expansion.

Using $\bar{\bm{u}}_{2}=\bm{x}_{s_{2}}+\mathcal{O}(h^{3})$ and $\lambda_{s_{2}}-\lambda_{s}=r_{2}h=\frac{2}{3}h$, we find that

Comparing with the Taylor expansion in Eq. (B.4) with $n=2$:

we need to check the following conditions:

The first two conditions are clear. The last condition follows from

Therefore, $\bar{\bm{x}}_{t}=\bm{x}_{t}+\mathcal{O}(h^{4})$. ∎

B.6 Connections to Explicit Exponential Runge-Kutta (expRK) Methods

Assume we have an ODE with the following form:

Appendix C Algorithms of DPM-Solvers

We firstly list the detailed DPM-Solver-1, 2, 3 in Algorithms 3, 4, 5. Note that DPM-Solver-2 is the general case with $r_{1}\in(0,1)$, and we usually set $r_{1}=0.5$ for DPM-Solver-2, as in Sec. 3.

Then we list the adaptive step size algorithms, named as DPM-Solver-12 (combining 1 and 2; Algorithm 6) and DPM-Solver-23 (combining 2 and 3; Algorithm 7). We follow to let the absolute tolerance $\epsilon_{\text{atol}}=\frac{\bm{x}_{\text{max}}-\bm{x}_{\text{min}}}{256}$ for image data, which is $0.0078$ for VP type DPMs. We can tune the relative tolerance $\epsilon_{\text{rtol}}$ to balance the accuracy and NFE, and we find that $\epsilon_{\text{rtol}}=0.05$ is good enough and can converge quickly.

In practice, the inputs of the adaptive step size solvers are batch data. We simply choose $E_{2}$ and $E_{3}$ as the maximum value of all the batch data. Besides, we implement the comparison $s>\epsilon$ by $|s-\epsilon|>10^{-5}$ to avoid numerical issues.

Appendix D Implementation Details of DPM-Solver

Theoretically, we need to solve diffusion ODEs from time $T$ to time to generate samples. Practically, the training and evaluation for the noise prediction model $\bm{\epsilon}_{\theta}(\bm{x}_{t},t)$ usually start from time $T$ to time $\epsilon$ to avoid numerical issues for $t$ near to , where $\epsilon>0$ is a hyperparameter .

In contrast to the sampling methods based on diffusion SDEs , We don’t add the “denoising” trick at the final step at time $\epsilon$ (which is to set the noise variance to zero), and we just solve diffusion ODEs from $T$ to $\epsilon$ by DPM-Solver, since we find it performs well enough.

For discrete-time DPMs, we firstly convert the model to continuous time (see Appendix D.2), and then solver it from time $T$ to time $t$.

D.2 Sampling from Discrete-Time DPMs

In this section, we discuss the more general case for discrete-time DPMs, in which we consider the 1000-step DPMs and the 4000-step DPMs , and we also consider the end time $\epsilon$ for sampling.

Type-1. Scale the discrete time steps $[t_{1},t_{N}]=[\frac{T}{N},T]$ to the continuous time range $[\frac{T}{N},T]$, and let $\bm{\epsilon}_{\theta}(\cdot,t)=\bm{\epsilon}_{\theta}(\cdot,\frac{T}{N})$ for $t\in[\epsilon,\frac{T}{N}]$. In this case, we can define the continuous-time noise prediction model by

where the continuous time $t\in[\epsilon,\frac{T}{N}]$ maps to the discrete input , and the continuous time $T$ maps to the discrete input $\frac{1000(N-1)}{N}$.

Type-2. Scale the discrete time steps $[t_{1},t_{N}]=[\frac{T}{N},T]$ to the continuous time range $[0,T]$. In this case, we can define the continuous-time noise prediction model by

where the continuous time maps to the discrete input , and the continuous time $T$ maps to the discrete input $\frac{1000(N-1)}{N}$.

In practice, we have $T=1$, and the smallest discrete time $t_{1}=10^{-3}$. For fixed $K$ number of function evaluations, we empirically find that for small $K$, the Type-1 with $\epsilon=10^{-3}$ may have better sample quality, and for large $K$, the Type-2 with $\epsilon=10^{-4}$ may have better sample quality. We refer to Appendix E for detailed results.

D.3 DPM-Solver in 20 Function Evaluations

Given a fixed budget $K\leq 20$ of the number of function evaluations, we uniformly divide the interval $[\lambda_{T},\lambda_{\epsilon}]$ into $M=(\lfloor K/3\rfloor+1)$ segments, and take $M$ steps to generate samples. The $M$ steps are dependent on the remainder $R$ of $K$ mod $3$ to make sure the total number of function evaluations is exactly $K$.

If $R=0$, we firstly take $M-2$ steps of DPM-Solver-3, and then take $1$ step of DPM-Solver-2 and 1 step of DPM-Solver-1. The total number of function evaluations is $3\cdot(\frac{K}{3}-1)+2+1=K$.

If $R=1$, we firstly take $M-1$ steps of DPM-Solver-3 and then take $1$ step of DPM-Solver-1. The total number of function evaluations is $3\cdot(\frac{K-1}{3})+1=K$.

If $R=2$, we firstly take $M-1$ steps of DPM-Solver-3 and then take $1$ step of DPM-Solver-2. The total number of function evaluations is $3\cdot(\frac{K-2}{3})+2=K$.

We empirically find that this design of time steps can greatly improve the generation quality, and DPM-Solver can generate comparable samples in 10 steps and high-quality samples in 20 steps.

The costs of computing $t_{\lambda}(\cdot)$ is negligible, because for the noise schedules of $\alpha_{t}$ and $\sigma_{t}$ used in previous DPMs (“linear” and “cosine”) , both $\lambda(t)$ and its inverse function $t_{\lambda}(\cdot)$ have analytic formulations. We mainly consider the variance preserving type here, since it is the most widely-used type. The functions of other types (variance exploding and sub-variance preserving type) can be similarly derived.

where $\beta_{0}=0.1$ and $\beta_{1}=20$, following . As $\sigma_{t}=\sqrt{1-\alpha_{t}^{2}}$, we can compute $\lambda_{t}$ analytically. Moreover, the inverse function is

To reduce the influence of numerical issues, we can compute $t_{\lambda}$ by the following equivalent formulation:

And we solve diffusion ODEs between $[\epsilon,T]$, where $T=1$.

where $s=0.008$, following . As clipped the derivatives to ensure the numerical stability, we also clip the maximum time $T=0.9946$. As $\sigma_{t}=\sqrt{1-\alpha_{t}^{2}}$, we can compute $\lambda_{t}$ analytically. Moreover, given a fixed $\lambda$, let

which computes the corresponding $\log\alpha$ for $\lambda$. Then the inverse function is

And we solve diffusion ODEs between $[\epsilon,T]$, where $T=0.9946$.

D.5 Conditional Sampling by DPM-Solver

DPM-Solver can also be used for conditional sampling, with a simple modification. The conditional generation needs to sample from the conditional diffusion ODE which includes the conditional noise prediction model. We follow the classifier guidance method to define the conditional noise prediction model as $\bm{\epsilon}_{\theta}(\bm{x}_{t},t,y)\coloneqq\bm{\epsilon}_{\theta}(\bm{x}_{t},t)-s\cdot\sigma_{t}\nabla_{\bm{x}}\log p_{t}(y|\bm{x}_{t};\theta)$, where $p_{t}(y|\bm{x}_{t};\theta)$ is a pre-trained classifier and $s$ is the classifier guidance scale (default is 1.0). Thus, we can use DPM-Solver to solve this diffusion ODE for fast conditional sampling, as shown in Fig. 1.

D.6 Numerical Stability

As we need to compute $e^{h_{i}}-1$ in the algorithm of DPM-Solver, we follow to use expm1($h_{i}$) instead of exp($h_{i}$)-1 to improve numerical stability.

Appendix E Experiment Details

For all experiments, we evaluate DPM-Solver on NVIDIA A40 GPUs. However, the computation resource can be other types of GPU, such as NVIDIA GeForce RTX 2080Ti, because we can tune the batch size for sampling.

Alternatively, the diffusion ODE can be reparameterized to the $\lambda$ domain. In this section, we propose the formulation of diffusion ODEs w.r.t. $\lambda$ for VP type, and other types can be similarly derived.

For a given $\lambda$, denote $\hat{\alpha}_{\lambda}\coloneqq\alpha_{t(\lambda)}$, $\hat{\sigma}_{\lambda}\coloneqq\sigma_{t(\lambda)}$. As $\hat{\alpha}_{\lambda}^{2}+\hat{\sigma}_{\lambda}^{2}=1$, we can prove that $\frac{\differential\lambda}{\differential\hat{\alpha}_{\lambda}}=\frac{1}{\hat{\alpha}_{\lambda}\hat{\sigma}^{2}_{\lambda}}$, so $\frac{\differential\log\hat{\alpha}_{\lambda}}{\differential\lambda}=\hat{\sigma}^{2}_{\lambda}$. Applying change-of-variable to Eq. (2.7), we have

The ODE Eq. (E.1) can be also solved directly by RK methods, and we use such formulation for the experiments of RK2 ($\lambda$) and RK3 ($\lambda$) in Table 1.

E.2 Code Implementation

We implement our code with both JAX (for continuous-time DPMs) and PyTorch (for discrete-time DPMs), and our code is released at https://github.com/LuChengTHU/dpm-solver.

E.3 Sample Quality Comparison with Continuous-Time Sampling Methods

Table 3 shows the detailed FID results, which is corresponding to Fig. 2(a). We use the official code and checkpoint in , the code license is Apache License 2.0. We use their released “checkpoint_8” of the “VP deep” type. We compare methods for $\epsilon=10^{-3}$ and $\epsilon=10^{-4}$. We find that the sampling methods based on diffusion SDEs can achieve better sample quality with $\epsilon=10^{-3}$; and that the sampling methods based on diffusion ODEs can achieve better sample quality with $\epsilon=10^{-4}$. For DPM-Solver, we find that DPM-Solver with less than 15 NFE can achieve better FID with $\epsilon=10^{-3}$ than $\epsilon=10^{-4}$, while DPM-Solver with more than 15 NFE can achieve better FID with $\epsilon=10^{-4}$ than $\epsilon=10^{-3}$.

For the diffusion SDEs with Euler discretization, we use the PC sampler in with “euler_maruyama” predictor and no corrector, which uses uniform time steps between $T$ and $\epsilon$. We add the “denoise” trick at the final step, which can greatly improve the FID score for $\epsilon=10^{-3}$.

For the diffusion SDEs with Improved Euler discretization , we follow the results in their original paper, which only includes the results with $\epsilon=10^{-3}$. The corresponding relative tolerance $\epsilon_{rel}$ are $0.50$, $0.10$ and $0.05$, respectively.

For the diffusion ODEs with RK45 Solver, we use the code in , and tune the atol and rtol of the solver. For the NFE from small to large, we use the same atol = rtol = $0.1$, $0.01$, $0.001$ for the results of $\epsilon=10^{-3}$, and the same atol = rtol = $0.1$, $0.05$, $0.02$, $0.01$, $0.001$ for the results of $\epsilon=10^{-4}$, respectively.

For the diffusion ODEs with DPM-Solver, we use the method in Appendix D.3 for NFE $\leq 20$, and the adaptive step size solver in Appendix C. For $\epsilon=10^{-3}$, we use DPM-Solver-12 with relative tolerance $\epsilon_{\text{rtol}}=0.05$. For $\epsilon=10^{-4}$, we use DPM-Solver-23 with relative tolerance $\epsilon_{\text{rtol}}=0.05$.

E.4 Sample Quality Comparison with RK Methods

Table 1 shows the different performance of RK methods and DPM-Solver-2 and 3. We list the detailed settings in this section.

We use $\bm{F}(\bm{x}_{t},t)=\bm{h}_{\theta}(\bm{x}_{t},t)$ in Eq. (2.7) for the results with RK2 ($t$) and RK3 ($t$), and $\bm{F}(\hat{\bm{x}}_{\lambda},\lambda)=\hat{\bm{h}}_{\theta}(\hat{\bm{x}}_{\lambda},\lambda)$ in Eq. (E.1) for the results with RK2 ($\lambda$) and RK3 ($\lambda$). For all experiments, we use the uniform step size w.r.t. $t$ or $\lambda$.

E.5 Sample Quality Comparison with Discrete-Time Sampling Methods

We compare DPM-Solver with other discrete-time sampling methods for DPMs, as shown in Table 4 and Table 5. We use the code in for sampling with DDPM and DDIM, and the code license is MIT License. We use the code in for sampling with Analytic-DDPM and Analytic-DDIM, whose license is unknown. We directly follow the best results in the original paper of GGDM .

For the CIFAR-10 experiments, we use the pretrained checkpoint by , which is also provided in the released code in . We use quadratic time steps for DDPM and DDIM, which empirically has better FID performance than the uniform time steps . We use the uniform time steps for Analytic-DDPM and Analytic-DDIM. For DPM-Solver, we use both Type-1 discrete and Type-2 discrete methods to convert the discrete-time model to the continuous-time model. We use the method in Appendix D.3 for NFE $\leq 20$, and the adaptive step size solver in Appendix C for NFE $>20$. For all the experiments, we use DPM-Solver-12 with relative tolerance $\epsilon_{\text{rtol}}=0.05$.

For the CelebA 64x64 experiments, we use the pretrained checkpoint by . We use quadratic time steps for DDPM and DDIM, which empirically has better FID performance than the uniform time steps . We use the uniform time steps for Analytic-DDPM and Analytic-DDIM. For DPM-Solver, we use both Type-1 discrete and Type-2 discrete methods to convert the discrete-time model to the continuous-time model. We use the method in Appendix D.3 for NFE $\leq 20$, and the adaptive step size solver in Appendix C for NFE $>20$. For all the experiments, we use DPM-Solver-12 with relative tolerance $\epsilon_{\text{rtol}}=0.05$. Note that our best FID results on CelebA 64x64 is even better than the 1000-step DDPM (and all the other methods).

For the ImageNet 64x64 experiments, we use the pretrained checkpoint by , and the code license is MIT License. We use the uniform time steps for DDPM and DDIM, following . We use the uniform time steps for Analytic-DDPM and Analytic-DDIM. For DPM-Solver, we use both Type-1 discrete and Type-2 discrete methods to convert the discrete-time model to the continuous-time model. We use the method in Appendix D.3 for NFE $\leq 20$, and the adaptive step size solver in Appendix C for NFE $>20$. For all the experiments, we use DPM-Solver-23 with relative tolerance $\epsilon_{\text{rtol}}=0.05$. Note that the ImageNet dataset includes real human photos and it may have privacy issues, as discussed in .

For the ImageNet 128x128 experiments, we use classifier guidance for sampling with the pretrained checkpoints (for both the diffusion model and the classifier model) by , and the code license is MIT License. We use the uniform time steps for DDPM and DDIM, following . For DPM-Solver, we only use Type-1 discrete method to convert the discrete-time model to the continuous-time model. We use the method in Appendix D.3 for NFE $\leq 20$, and the adaptive step size solver DPM-Solver-12 with relative tolerance $\epsilon_{\text{rtol}}=0.05$ (detailed in Appendix C) for NFE $>20$. For all the experiments, we set the classifier guidance scale $s=1.25$, which is the best setting for DDIM in (we refer to their Table 14 for details).

For the LSUN bedroom 256x256 experiments, we use the unconditional pretrained checkpoint by , and the code license is MIT License. We use the uniform time steps for DDPM and DDIM, following . For DPM-Solver, we only use Type-1 discrete method to convert the discrete-time model to the continuous-time model. We use the method in Appendix D.3 for DPM-Solver.

E.6 Comparing Different Orders of DPM-Solver

We also compare the sample quality of the different orders of DPM-Solver, as shown in Table 6. We use DPM-Solver-1,2,3 with uniform time steps w.r.t. $\lambda$, and the fast version in Appendix D.3 for NFE less than 20, and we name it as DPM-Solver-fast. For the discrete-time models, we only compare the Type-2 discrete method, and the results of Type-1 are similar.

As the actual NFE of DPM-Solver-2 is $2\times\lfloor\text{NFE}/2\rfloor$ and the actual NFE of DPM-Solver-3 is $3\times\lfloor\text{NFE}/3\rfloor$, which may be smaller than NFE, we use the notation † to note that the actual NFE is less than the given NFE. We find that for NFE less than 20, the proposed fast version (DPM-Solver-fast) is usually better than the single order method, and for larger NFE, DPM-Solver-3 is better than DPM-Solver-2, and DPM-Solver-2 is better than DPM-Solver-1, which matches our proposed convergence rate analysis.

E.7 Runtime Comparison between DPM-Solver and DDIM

Theoretically, for the same NFE, the runtime of DPM-Solver and DDIM are almost the same (linear to NFE) because the main computation costs are the serial evaluations of the large neural network $\bm{\epsilon}_{\theta}$ and the other coefficients are analytically computed with ignorable costs.

Table 7 shows the runtime of DPM-Solver and DDIM on a single NVIDIA A40, varying different datasets and NFE. We use torch.cuda.Event and torch.cuda.synchronize for accurately computing the runtime. We use the discrete-time pretrained diffusion models for each dataset. We evaluate the runtime for 8 batches and computes the mean and std of the runtime. We use 64 batch size for LSUN bedroom 256x256 due to the GPU memory limitation, and 128 batch size for other datasets.

For DDIM, we use the official implementationhttps://github.com/ermongroup/ddim. We find that our implementation of DPM-Solver reduces some repetitive computation of the coefficients, so under the same NFE, DPM-Solver is slightly faster than DDIM of their implementation. Nevertheless, the runtime evaluation results show that the runtime of DPM-Solver and DDIM are almost the same for the same NFE, and the runtime is approximately linear to the NFE. Therefore, the speedup for the NFE is almost the actual speedup of the runtime, so the proposed DPM-Solver can greatly speedup the sampling of DPMs.

E.8 Conditional Sampling on ImageNet 256x256

For the conditional sampling in Fig. 1, we use the pretrained checkpoint in with classifier guidance (ADM-G), and the classifier scale is $1.0$. The code license is MIT License. We use uniform time step for DDIM, and the fast version for DPM-Solver in Appendix D.3 (DPM-Solver-fast) with 10, 15, 20 and 100 steps.

Fig. 3 shows the conditional sample results by DDIM and DPM-Solver. We find that DPM-Solver with 15 NFE can generate comparable samples with DDIM with 100 NFE.

E.9 Additional Samples

Additional sampling results on CIFAR-10, CelebA 64x64, ImageNet 64x64, LSUN bedroom 256x256 , ImageNet 256x256 are reported in Figs. 4-8.