diff --git a/.github/workflows/bump.yml b/.github/workflows/bump.yml
index d8639b93..1995fa31 100644
--- a/.github/workflows/bump.yml
+++ b/.github/workflows/bump.yml
@@ -27,6 +27,9 @@ jobs:
     - name: Bump version
       run: git branch --show-current | sed 's|release/||' | xargs poetry version | { printf '::set-output name=PR_TITLE::'; cat; }
       id: bump
+    
+    - name: Bump pfhedge.__version__
+      run:  git branch --show-current | sed 's|release/||' | xargs -I {} echo '__version__ = "{}"' > pfhedge/version.py
 
     - name: Create pull request
       uses: peter-evans/create-pull-request@v3
diff --git a/.github/workflows/publish.yml b/.github/workflows/publish.yml
index b6a58207..298fdc55 100644
--- a/.github/workflows/publish.yml
+++ b/.github/workflows/publish.yml
@@ -1,6 +1,10 @@
 name: Publish
 
 on:
+  push:
+    branches: main
+  release:
+    types: [published]
   workflow_dispatch:
 
 jobs:
diff --git a/docs/source/nn.functional.rst b/docs/source/nn.functional.rst
index f507ed5d..e842cc28 100644
--- a/docs/source/nn.functional.rst
+++ b/docs/source/nn.functional.rst
@@ -19,6 +19,16 @@ Payoff Functions
     european_binary_payoff
     european_forward_start_payoff
 
+Nonlinear activation functions
+------------------------------
+
+.. autosummary::
+    :nosignatures:
+    :toctree: generated
+
+    leaky_clamp
+    clamp
+
 Criterion Functions
 -------------------
 
@@ -32,21 +42,48 @@ Criterion Functions
     expected_shortfall
     value_at_risk
 
+Black-Scholes formulas
+----------------------
+
+.. autosummary::
+    :nosignatures:
+    :toctree: generated
+
+    bs_european_price
+    bs_european_delta
+    bs_european_gamma
+    bs_european_vega
+    bs_european_theta
+    bs_american_binary_price
+    bs_american_binary_delta
+    bs_american_binary_gamma
+    bs_american_binary_vega
+    bs_american_binary_theta
+    bs_european_binary_price
+    bs_european_binary_delta
+    bs_european_binary_gamma
+    bs_european_binary_vega
+    bs_european_binary_theta
+    bs_lookback_price
+    bs_lookback_delta
+    bs_lookback_gamma
+    bs_lookback_vega
+    bs_lookback_theta
+
 Other Functions
------------------
+---------------
 
 .. autosummary::
     :nosignatures:
     :toctree: generated
 
-    leaky_clamp
-    clamp
-    topp
-    realized_variance
-    realized_volatility
-    terminal_value
-    ncdf
-    npdf
+    bilerp
     d1
     d2
+    ncdf
+    npdf
+    realized_variance
+    realized_volatility
     svi_variance
+    terminal_value
+    topp
diff --git a/docs/source/notes/adding_clause.rst b/docs/source/notes/adding_clause.rst
new file mode 100644
index 00000000..dbee6b56
--- /dev/null
+++ b/docs/source/notes/adding_clause.rst
@@ -0,0 +1,69 @@
+Adding clause to derivative
+===========================
+
+One can customize a derivative by registering additional clauses;
+see :func:`pfhedge.instruments.BaseDerivative.add_clause`.
+
+Let us see how to add clauses to derivatives by taking European capped call option as an example.
+This option is a variant of a European option, which is given by
+
+.. code:: python
+
+    >>> from pfhedge.instruments import BrownianStock
+    >>> from pfhedge.instruments import EuropeanOption
+    ...
+    >>> strike = 1.0
+    >>> maturity = 1.0
+    >>> stock = BrownianStock()
+    >>> european = EuropeanOption(stock, strike=strike, maturity=maturity)
+
+The capped call associates ''a barrier clause'';
+if the underlier's spot reaches the barrier price :math:`B`,
+being greater than the strike :math:`K` and the spot at inception,
+the option immediately expires and pays off its intrinsic value at that moment :math:`B - K`.
+
+This clause can be registered to a derivative as follows:
+
+.. code:: python
+
+    >>> def cap_clause(derivative, payoff):
+    ...     barrier = 1.4
+    ...     max_spot = derivative.ul().spot.max(-1).values
+    ...     capped_payoff = torch.full_like(payoff, barrier - strike)
+    ...     return torch.where(max_spot < barrier, payoff, capped_payoff)
+    ...
+    >>> capped_european = EuropeanOption(stock, strike=strike, maturity=maturity)
+    >>> capped_european.add_clause("cap_clause", cap_clause)
+
+The method ``add_clause`` adds the clause and its name to the derivative.
+Here the function ``cap_caluse`` represents the clause to modify the payoff depending on the state of the derivative.
+The clause function should have the signature ``clause(derivative, payoff) -> modified payoff``.
+
+The payoff would be capped as intended:
+
+.. code:: python
+
+    >>> n_paths = 100000
+    >>> capped_european.simulate(n_paths=n_paths)
+    >>> european.payoff().max()
+    >>> # 1.2...
+    >>> capped_european.payoff().max()
+    >>> # 0.4
+
+The price of the capped European call option can be evaluated by using a European option as a control variates.
+
+.. code:: python
+
+    >>> from math import sqrt
+    ...
+    >>> payoff_european = european.payoff()
+    >>> payoff_capped_european = capped_european.payoff()
+    >>> bs_price = BlackScholes(european).price(0.0, european.maturity, stock.sigma).item()
+    >>> price = bs_price + (payoff_capped_european - payoff_european).mean().item()
+    >>> error = (payoff_capped_european - payoff_european).std().item() / sqrt(n_paths)
+    >>> bs_price
+    >>> # 0.07967...
+    >>> price
+    >>> # 0.07903...
+    >>> error
+    >>> # 0.00012...
diff --git a/examples/example_adding_clause.py b/examples/example_adding_clause.py
new file mode 100644
index 00000000..c5a6ff5a
--- /dev/null
+++ b/examples/example_adding_clause.py
@@ -0,0 +1,57 @@
+import sys
+
+sys.path.append("..")
+
+from math import sqrt
+
+import torch
+
+from pfhedge.instruments import BrownianStock
+from pfhedge.instruments import EuropeanOption
+from pfhedge.nn import BlackScholes
+
+
+def main():
+    torch.manual_seed(42)
+
+    strike = 1.0
+    maturity = 1.0
+    stock = BrownianStock()
+    european = EuropeanOption(stock, strike=strike, maturity=maturity)
+
+    def cap_clause(derivative, payoff):
+        barrier = 1.4
+        max_spot = derivative.ul().spot.max(-1).values
+        capped_payoff = torch.full_like(payoff, barrier - strike)
+        return torch.where(max_spot < barrier, payoff, capped_payoff)
+
+    capped_european = EuropeanOption(stock, strike=strike, maturity=maturity)
+    capped_european.add_clause("cap_clause", cap_clause)
+
+    n_paths = 100000
+    capped_european.simulate(n_paths=n_paths)
+
+    payoff_european = european.payoff()
+    payoff_capped_european = capped_european.payoff()
+
+    max_spot = payoff_european.max().item()
+    capped_max_spot = payoff_capped_european.max().item()
+    print("Max payoff of vanilla European:", max_spot)
+    print("Max payoff of capped  European:", capped_max_spot)
+
+    # Price using control variates
+    bs_price = BlackScholes(european).price(0.0, european.maturity, stock.sigma).item()
+    value0 = payoff_capped_european.mean().item()
+    value1 = bs_price + (payoff_capped_european - payoff_european).mean().item()
+    error0 = payoff_capped_european.std().item() / sqrt(n_paths)
+    error1 = (payoff_capped_european - payoff_european).std().item() / sqrt(n_paths)
+
+    print("BS price of vanilla European:", bs_price)
+    print("Price of capped European without control variates:", value0)
+    print("Price of capped European with    control variates:", value1)
+    print("Error of capped European without control variates:", error0)
+    print("Error of capped European with    control variates:", error1)
+
+
+if __name__ == "__main__":
+    main()
diff --git a/pfhedge/__init__.py b/pfhedge/__init__.py
index 29da3648..a267f9be 100644
--- a/pfhedge/__init__.py
+++ b/pfhedge/__init__.py
@@ -1,3 +1,5 @@
 # Users can import Hedger either as `from pfhedge import Hedger` or
 # `from pfhedge.nn import Hedger`
 from pfhedge.nn import Hedger
+
+from .version import __version__
diff --git a/pfhedge/_utils/bisect.py b/pfhedge/_utils/bisect.py
index f2010810..44da9fcf 100644
--- a/pfhedge/_utils/bisect.py
+++ b/pfhedge/_utils/bisect.py
@@ -1,3 +1,4 @@
+from typing import Any
 from typing import Callable
 from typing import Union
 
@@ -86,7 +87,7 @@ def find_implied_volatility(
     upper: float = 1.000,
     precision: float = 1e-6,
     max_iter: int = 100,
-    **params,
+    **params: Any,
 ) -> Tensor:
     """Find implied volatility by binary search.
 
diff --git a/pfhedge/_utils/operations.py b/pfhedge/_utils/operations.py
index bce4952f..bb67663d 100644
--- a/pfhedge/_utils/operations.py
+++ b/pfhedge/_utils/operations.py
@@ -1,3 +1,4 @@
+from typing import Any
 from typing import Callable
 
 import torch
@@ -5,7 +6,7 @@
 
 
 def ensemble_mean(
-    function: Callable[..., Tensor], n_times: int = 1, *args, **kwargs
+    function: Callable[..., Tensor], n_times: int = 1, *args: Any, **kwargs: Any
 ) -> Tensor:
     """Compute ensemble mean from function.
 
diff --git a/pfhedge/_utils/parse.py b/pfhedge/_utils/parse.py
index 1ed5c3b7..987e24bc 100644
--- a/pfhedge/_utils/parse.py
+++ b/pfhedge/_utils/parse.py
@@ -1,4 +1,5 @@
 from numbers import Real
+from typing import Any
 from typing import Optional
 from typing import Union
 
@@ -16,7 +17,7 @@ def parse_spot(
     strike: Optional[Tensor] = None,
     moneyness: Optional[Tensor] = None,
     log_moneyness: Optional[Tensor] = None,
-    **kwargs
+    **kwargs: Any
 ) -> Tensor:
     spot = _as_optional_tensor(spot)
     strike = _as_optional_tensor(strike)
@@ -34,7 +35,10 @@ def parse_spot(
 
 
 def parse_volatility(
-    *, volatility: Optional[Tensor] = None, variance: Optional[Tensor] = None, **kwargs
+    *,
+    volatility: Optional[Tensor] = None,
+    variance: Optional[Tensor] = None,
+    **kwargs: Any
 ) -> Tensor:
     if volatility is not None:
         return volatility
@@ -45,7 +49,7 @@ def parse_volatility(
 
 
 def parse_time_to_maturity(
-    *, time_to_maturity: Optional[Tensor] = None, **kwargs
+    *, time_to_maturity: Optional[Tensor] = None, **kwargs: Any
 ) -> Tensor:
     if time_to_maturity is not None:
         return time_to_maturity
diff --git a/pfhedge/_utils/testing.py b/pfhedge/_utils/testing.py
index 752ec6d1..798e6c72 100644
--- a/pfhedge/_utils/testing.py
+++ b/pfhedge/_utils/testing.py
@@ -1,3 +1,4 @@
+from typing import Any
 from typing import Callable
 
 import torch
@@ -62,7 +63,7 @@ def assert_convex(
 
 
 def assert_cash_invariant(
-    fn: Callable[[Tensor], Tensor], x: Tensor, c: float, **kwargs
+    fn: Callable[[Tensor], Tensor], x: Tensor, c: float, **kwargs: Any
 ) -> None:
     """Assert cash invariance.
 
@@ -78,7 +79,7 @@ def assert_cash_invariant(
 
 
 def assert_cash_equivalent(
-    fn: Callable[[Tensor], Tensor], x: Tensor, c: float, **kwargs
+    fn: Callable[[Tensor], Tensor], x: Tensor, c: float, **kwargs: Any
 ) -> None:
     """Assert ``c`` is the cash equivalent of ``x``.
     ``fn(x) = fn(torch.full_like(x, c))``
diff --git a/pfhedge/autogreek.py b/pfhedge/autogreek.py
index 2e0d5d88..62366854 100644
--- a/pfhedge/autogreek.py
+++ b/pfhedge/autogreek.py
@@ -1,4 +1,5 @@
 from inspect import signature
+from typing import Any
 from typing import Callable
 
 import torch
@@ -10,7 +11,7 @@
 
 
 def delta(
-    pricer: Callable[..., Tensor], *, create_graph: bool = False, **params
+    pricer: Callable[..., Tensor], *, create_graph: bool = False, **params: Any
 ) -> Tensor:
     """Computes and returns delta of a derivative using automatic differentiation.
 
@@ -42,12 +43,14 @@ def delta(
         >>> import pfhedge.autogreek as autogreek
         >>> from pfhedge.nn import BSEuropeanOption
         >>>
+        >>> # TODO(simaki): Rewrite using functional
         >>> pricer = BSEuropeanOption().price
         >>> autogreek.delta(
         ...     pricer,
         ...     log_moneyness=torch.zeros(3),
         ...     time_to_maturity=torch.ones(3),
         ...     volatility=torch.tensor([0.18, 0.20, 0.22]),
+        ...     strike=1.0,
         ... )
         tensor([0.5359, 0.5398, 0.5438])
 
@@ -81,11 +84,6 @@ def delta(
         >>> autogreek.delta(pricer, spot=torch.tensor(1.0))
         tensor(0.5...)
     """
-    if params.get("strike") is None and params.get("spot") is None:
-        # Since delta does not depend on strike,
-        # assign an arbitrary value (1.0) to strike if not given.
-        params["strike"] = torch.tensor(1.0)
-
     spot = parse_spot(**params).requires_grad_()
     params["spot"] = spot
     if "strike" in params:
@@ -108,7 +106,7 @@ def delta(
 
 
 def gamma(
-    pricer: Callable[..., Tensor], *, create_graph: bool = False, **params
+    pricer: Callable[..., Tensor], *, create_graph: bool = False, **params: Any
 ) -> Tensor:
     """Computes and returns gamma of a derivative.
 
@@ -134,7 +132,6 @@ def gamma(
         torch.Tensor
 
     Examples:
-
         Gamma of a European option can be evaluated as follows.
 
         >>> import pfhedge.autogreek as autogreek
@@ -175,7 +172,7 @@ def gamma(
 
 
 def vega(
-    pricer: Callable[..., Tensor], *, create_graph: bool = False, **params
+    pricer: Callable[..., Tensor], *, create_graph: bool = False, **params: Any
 ) -> Tensor:
     """Computes and returns vega of a derivative using automatic differentiation.
 
@@ -200,7 +197,6 @@ def vega(
         torch.Tensor
 
     Examples:
-
         Vega of a European option can be evaluated as follows.
 
         >>> import pfhedge.autogreek as autogreek
@@ -235,7 +231,7 @@ def vega(
 
 
 def theta(
-    pricer: Callable[..., Tensor], *, create_graph: bool = False, **params
+    pricer: Callable[..., Tensor], *, create_graph: bool = False, **params: Any
 ) -> Tensor:
     """Computes and returns theta of a derivative using automatic differentiation.
 
@@ -259,7 +255,6 @@ def theta(
         torch.Tensor
 
     Examples:
-
         Theta of a European option can be evaluated as follows.
 
         >>> import pfhedge.autogreek as autogreek
diff --git a/pfhedge/features/_getter.py b/pfhedge/features/_getter.py
index ff5f0884..643c2258 100644
--- a/pfhedge/features/_getter.py
+++ b/pfhedge/features/_getter.py
@@ -1,4 +1,5 @@
 from collections import OrderedDict
+from typing import Any
 from typing import Dict
 from typing import Iterator
 from typing import Tuple
@@ -74,7 +75,7 @@ def get_class(self, name: str) -> Type[Feature]:
             )
         return self._features[name]
 
-    def get_instance(self, name: str, *args, **kwargs) -> Feature:
+    def get_instance(self, name: str, **kwargs: Any) -> Feature:
         """Returns the feature with the given name.
 
         Parameters:
@@ -83,22 +84,21 @@ def get_instance(self, name: str, *args, **kwargs) -> Feature:
         Returns:
             Feature: feature.
         """
-        return self.get_class(name)(*args, **kwargs)  # type: ignore
+        return self.get_class(name)(**kwargs)  # type: ignore
 
 
-def get_feature(feature: Union[str, Feature], *args, **kwargs) -> Feature:
+def get_feature(feature: Union[str, Feature], **kwargs: Any) -> Feature:
     """Get feature from name.
 
     Args:
         name (str): Name of feature.
-        *args: Arguments to pass to feature constructor.
-        *kwargs: Keyword arguments to pass to feature constructor.
+        **kwargs: Keyword arguments to pass to feature constructor.
 
     Returns:
         Feature
     """
     if isinstance(feature, str):
-        feature = FeatureFactory().get_instance(feature, *args, **kwargs)
+        feature = FeatureFactory().get_instance(feature, **kwargs)
     elif not isinstance(feature, Feature):
         raise TypeError(f"{feature} is not an instance of Feature.")
     return feature
diff --git a/pfhedge/features/container.py b/pfhedge/features/container.py
index df38d6a5..cff67de7 100644
--- a/pfhedge/features/container.py
+++ b/pfhedge/features/container.py
@@ -14,6 +14,7 @@
 from ._getter import get_feature
 
 T = TypeVar("T", bound="FeatureList")
+TM = TypeVar("TM", bound="ModuleOutput")
 
 
 class FeatureList(Feature):
@@ -73,7 +74,6 @@ class ModuleOutput(Feature, Module):
         inputs (list[Feature]): The input features to the module.
 
     Examples:
-
         >>> from torch.nn import Linear
         >>> from pfhedge.instruments import BrownianStock
         >>> from pfhedge.instruments import EuropeanOption
@@ -121,7 +121,7 @@ def forward(self, input: Tensor) -> Tensor:
     def get(self, time_step: Optional[int]) -> Tensor:
         return self(self.inputs.get(time_step))
 
-    def of(self, derivative=None, hedger=None):
+    def of(self: TM, derivative=None, hedger=None) -> TM:
         self.inputs = self.inputs.of(derivative, hedger)
         return self
 
diff --git a/pfhedge/instruments/base.py b/pfhedge/instruments/base.py
index 584621c9..5157179a 100644
--- a/pfhedge/instruments/base.py
+++ b/pfhedge/instruments/base.py
@@ -1,5 +1,6 @@
 from abc import ABC
 from abc import abstractmethod
+from typing import Any
 from typing import List
 from typing import Optional
 from typing import TypeVar
@@ -42,7 +43,7 @@ def simulate(self, n_paths: int, time_horizon: float, **kwargs) -> None:
         """
 
     @abstractmethod
-    def to(self: T, *args, **kwargs) -> T:
+    def to(self: T, *args: Any, **kwargs: Any) -> T:
         """Moves and/or casts the buffers of the instrument.
 
         This can be called as
@@ -144,6 +145,18 @@ def half(self: T) -> T:
         """
         return self.to(torch.float16)
 
+    def float64(self: T) -> T:
+        """Alias for :meth:`double()`."""
+        return self.double()
+
+    def float32(self: T) -> T:
+        """Alias for :meth:`float()`."""
+        return self.float()
+
+    def float16(self: T) -> T:
+        """Alias for :meth:`half()`."""
+        return self.half()
+
     def bfloat16(self: T) -> T:
         """Casts all floating point parameters and buffers to
         ``torch.bfloat16`` datatype.
@@ -195,8 +208,8 @@ def _dinfo(self) -> List[str]:
 
 
 class Instrument(BaseInstrument):
-    def __init__(self, *args, **kwargs):
-        super().__init__(*args, **kwargs)
+    def __init__(self, *args, **kwargs) -> None:
+        super().__init__(*args, **kwargs)  # type: ignore
         raise DeprecationWarning(
             "Instrument is deprecated. Use BaseInstrument instead."
         )
diff --git a/pfhedge/instruments/derivative/american_binary.py b/pfhedge/instruments/derivative/american_binary.py
index 4ca3e912..039990e3 100644
--- a/pfhedge/instruments/derivative/american_binary.py
+++ b/pfhedge/instruments/derivative/american_binary.py
@@ -14,15 +14,7 @@
 
 
 class AmericanBinaryOption(BaseOption):
-    r"""American binary Option.
-
-    An American binary call option pays an unit amount of cash if and only if
-    the maximum of the underlying asset's price until maturity is equal or greater
-    than the strike price.
-
-    An American binary put option pays an unit amount of cash if and only if
-    the minimum of the underlying asset's price until maturity is equal or smaller
-    than the strike price.
+    r"""American binary option.
 
     The payoff of an American binary call option is given by:
 
@@ -31,10 +23,11 @@ class AmericanBinaryOption(BaseOption):
         \begin{cases}
             1 & (\mathrm{Max} \geq K) \\
             0 & (\text{otherwise})
-        \end{cases}
+        \end{cases} ,
 
-    Here, :math:`\mathrm{Max}` is the maximum of the underlying asset's price
-    until maturity and :math:`K` is the strike price (`strike`) of the option.
+    where
+    :math:`\mathrm{Max}` is the maximum of the underlier's spot price until maturity and
+    :math:`K` is the strike.
 
     The payoff of an American binary put option is given by:
 
@@ -43,12 +36,13 @@ class AmericanBinaryOption(BaseOption):
         \begin{cases}
             1 & (\mathrm{Min} \leq K) \\
             0 & (\text{otherwise})
-        \end{cases}
+        \end{cases} ,
 
-    Here, :math:`\mathrm{Min}` is the minimum of the underlying asset's price.
+    where
+    :math:`\mathrm{Min}` is the minimum of the underlier's spot price until maturity.
 
     .. seealso::
-        :func:`pfhedge.nn.functional.american_binary_payoff`: Payoff function.
+        - :func:`pfhedge.nn.functional.american_binary_payoff`
 
     Args:
         underlier (:class:`BasePrimary`): The underlying instrument of the option.
diff --git a/pfhedge/instruments/derivative/base.py b/pfhedge/instruments/derivative/base.py
index a716693e..41ceedbb 100644
--- a/pfhedge/instruments/derivative/base.py
+++ b/pfhedge/instruments/derivative/base.py
@@ -3,6 +3,7 @@
 from typing import Any
 from typing import Callable
 from typing import Dict
+from typing import Iterator
 from typing import Optional
 from typing import Tuple
 from typing import TypeVar
@@ -19,6 +20,7 @@
 from ..primary.base import BasePrimary
 
 T = TypeVar("T", bound="BaseDerivative")
+Clause = Callable[[T, Tensor], Tensor]
 
 
 class BaseDerivative(BaseInstrument):
@@ -80,7 +82,7 @@ def ul(self) -> BasePrimary:
         """Alias for ``self.underlier``."""
         return self.underlier
 
-    def to(self: T, *args, **kwargs) -> T:
+    def to(self: T, *args: Any, **kwargs: Any) -> T:
         self.underlier.to(*args, **kwargs)
         return self
 
@@ -116,7 +118,7 @@ def payoff(self) -> Tensor:
             torch.Tensor
         """
         payoff = self.payoff_fn()
-        for clause in self._clauses.values():
+        for clause in self.clauses():
             payoff = clause(self, payoff)
         return payoff
 
@@ -142,14 +144,13 @@ def delist(self: T) -> None:
         After this method self will be a private derivative.
         """
         self.pricer = None
+        self.cost = 0.0
 
     @property
     def is_listed(self) -> bool:
         return self.pricer is not None
 
-    def add_clause(
-        self, name: str, clause: Callable[["BaseDerivative", Tensor], Tensor]
-    ) -> None:
+    def add_clause(self, name: str, clause: Clause) -> None:
         """Adds a clause to the derivative.
 
         The clause will be called after :meth:`payoff_fn` method
@@ -173,8 +174,23 @@ def add_clause(
             raise KeyError('clause name can\'t contain ".", got: {}'.format(name))
         elif name == "":
             raise KeyError('clause name can\'t be empty string ""')
+
+        if not hasattr(self, "_clauses"):
+            raise AttributeError(
+                "cannot assign clause before BaseDerivative.__init__() call"
+            )
+
         self._clauses[name] = clause
 
+    def named_clauses(self) -> Iterator[Tuple[str, Clause]]:
+        if hasattr(self, "_clauses"):
+            for name, clause in self._clauses.items():
+                yield name, clause
+
+    def clauses(self) -> Iterator[Clause]:
+        for _, clause in self.named_clauses():
+            yield clause
+
     @property
     def spot(self) -> Tensor:
         """Returns ``self.pricer(self)`` if self is listed.
@@ -198,8 +214,8 @@ def __repr__(self) -> str:
 
 
 class Derivative(BaseDerivative):
-    def __init__(self, *args, **kwargs):
-        super().__init__(*args, **kwargs)
+    def __init__(self, *args, **kwargs) -> None:
+        super().__init__(*args, **kwargs)  # type: ignore
         raise DeprecationWarning(
             "Derivative is deprecated. Use BaseDerivative instead."
         )
diff --git a/pfhedge/instruments/derivative/cliquet.py b/pfhedge/instruments/derivative/cliquet.py
index d1e11115..b9a3406c 100644
--- a/pfhedge/instruments/derivative/cliquet.py
+++ b/pfhedge/instruments/derivative/cliquet.py
@@ -21,15 +21,18 @@ class EuropeanForwardStartOption(BaseDerivative):
         \mathrm{payoff} = \max(S_T / S_{T'} - K, 0) ,
 
     where
-    :math:`S_T` is the underlying asset's price at maturity,
-    :math:`S_{T'}` is the underlying asset's price at start, and
-    :math:`K` is the strike value of the option.
+    :math:`S_T` is the underlier's spot price at maturity,
+    :math:`S_{T'}` is the underlier's spot price at ``start``, and
+    :math:`K` is the strike.
 
     Note:
         If ``start`` is not divisible by the interval of the time step,
         it is rounded off so that the start time is the largest value
         that is divisible by the interval and less than ``start``.
 
+    .. seealso::
+        - :func:`pfhedge.nn.functional.european_forward_start_payoff`
+
     Args:
         underlier (:class:`BasePrimary`): The underlying instrument of the option.
         strike (float, default=1.0): The strike value of the option.
diff --git a/pfhedge/instruments/derivative/european.py b/pfhedge/instruments/derivative/european.py
index 8009b488..21f26342 100644
--- a/pfhedge/instruments/derivative/european.py
+++ b/pfhedge/instruments/derivative/european.py
@@ -16,25 +16,22 @@
 class EuropeanOption(BaseOption):
     r"""European option.
 
-    A European option provides its holder the right to buy (for call option)
-    or sell (for put option) an underlying asset with the strike price
-    on the date of maturity.
-
     The payoff of a European call option is given by:
 
     .. math::
-        \mathrm{payoff} = \max(S - K, 0)
+        \mathrm{payoff} = \max(S - K, 0) ,
 
-    Here, :math:`S` is the underlying asset's price at maturity and
-    :math:`K` is the strike price (`strike`) of the option.
+    where
+    :math:`S` is the underlier's spot price at maturity and
+    :math:`K` is the strike.
 
     The payoff of a European put option is given by:
 
     .. math::
-        \mathrm{payoff} = \max(K - S, 0)
+        \mathrm{payoff} = \max(K - S, 0) .
 
     .. seealso::
-        :func:`pfhedge.nn.functional.european_payoff`: Payoff function.
+        - :func:`pfhedge.nn.functional.european_payoff`
 
     Args:
         underlier (:class:`BasePrimary`): The underlying instrument of the option.
diff --git a/pfhedge/instruments/derivative/european_binary.py b/pfhedge/instruments/derivative/european_binary.py
index e252dfeb..86a72c24 100644
--- a/pfhedge/instruments/derivative/european_binary.py
+++ b/pfhedge/instruments/derivative/european_binary.py
@@ -16,41 +16,36 @@
 class EuropeanBinaryOption(BaseOption):
     r"""European binary option.
 
-    An American binary call option pays an unit amount of cash if and only if
-    the underlying asset's price at maturity is equal or greater than the strike price.
-
-    An American binary put option pays an unit amount of cash if and only if
-    the underlying asset's price at maturity is equal or smaller than the strike price.
-
-    The payoff of an American binary call option is given by:
+    The payoff of a European binary call option is given by
 
     .. math::
         \mathrm{payoff} =
         \begin{cases}
             1 & (S \geq K) \\
             0 & (\text{otherwise})
-        \end{cases}
+        \end{cases} ,
 
-    with :math:`S` being the underlying asset's price at maturity and
-    :math:`K` being the strike price (`strike`) of the option
+    where
+    :math:`S` is the underlier's spot price at maturity and
+    :math:`K` is the strike.
 
-    The payoff of an American binary put option is given by:
+    The payoff of a European binary put option is given by
 
     .. math::
         \mathrm{payoff} =
         \begin{cases}
             1 & (S \leq K) \\
             0 & (\text{otherwise})
-        \end{cases}
+        \end{cases} .
 
     .. seealso::
-        :func:`pfhedge.nn.functional.european_binary_payoff`: Payoff function.
+        - :func:`pfhedge.nn.functional.european_binary_payoff`
 
     Args:
         underlier (:class:`BasePrimary`): The underlying instrument of the option.
         call (bool, default=True): Specifies whether the option is call or put.
         strike (float, default=1): The strike price of the option.
-        maturity (float, default=20/250) The maturity of the option.
+        maturity (float, default=20/250): The maturity of the option.
 
     Attributes:
         dtype (torch.dtype): The dtype with which the simulated time-series are
diff --git a/pfhedge/instruments/derivative/lookback.py b/pfhedge/instruments/derivative/lookback.py
index 4d02e942..7a80025f 100644
--- a/pfhedge/instruments/derivative/lookback.py
+++ b/pfhedge/instruments/derivative/lookback.py
@@ -16,29 +16,25 @@
 class LookbackOption(BaseOption):
     r"""Lookback option with fixed strike.
 
-    A lookback call option provides its holder the right to buy an underlying with
-    the strike price and to sell with the highest price until the date of maturity.
-
-    A lookback put option provides its holder the right to sell an underlying with
-    the strike price and to buy with the lowest price until the date of maturity.
-
-    The payoff of a lookback call option is given by:
+    The payoff of a lookback call option is given by
 
     .. math::
-        \mathrm{payoff} = \max(\mathrm{Max} - K, 0)
+        \mathrm{payoff} = \max(\mathrm{Max} - K, 0) ,
 
-    Here, :math:`\mathrm{Max}` is the maximum of the underlying asset's price
-    until maturity and :math:`K` is the strike price (`strike`) of the option.
+    where
+    :math:`\mathrm{Max}` is the maximum of the underlier's spot price until maturity and
+    :math:`K` is the strike.
 
-    The payoff of a lookback put option is given by:
+    The payoff of a lookback put option is given by
 
     .. math::
-        \mathrm{payoff} = \max(K - \mathrm{Min}, 0)
+        \mathrm{payoff} = \max(K - \mathrm{Min}, 0) ,
 
-    Here, :math:`\mathrm{Min}` is the minimum of the underlying asset's price.
+    where
+    :math:`\mathrm{Max}` is the minimum of the underlier's spot price until maturity.
 
     .. seealso::
-        :func:`pfhedge.nn.functional.lookback_payoff`: Payoff function.
+        - :func:`pfhedge.nn.functional.lookback_payoff`
 
     Args:
         underlier (:class:`BasePrimary`): The underlying instrument of the option.
diff --git a/pfhedge/instruments/derivative/variance_swap.py b/pfhedge/instruments/derivative/variance_swap.py
index 37e8f0b0..0d61a4df 100644
--- a/pfhedge/instruments/derivative/variance_swap.py
+++ b/pfhedge/instruments/derivative/variance_swap.py
@@ -15,19 +15,21 @@
 class VarianceSwap(BaseDerivative):
     r"""Variance swap.
 
-    A variance swap pays cash in the amount of the realized variance
-    until the maturity and levies the cash of the strike variance.
-
     The payoff of a variance swap is given by
 
     .. math::
 
-        \mathrm{payoff} = \sigma^2 - K
+        \mathrm{payoff} = \sigma^2 - K ,
+
+    where
+    :math:`\sigma^2` is the realized variance of the underlying asset until maturity and
+    :math:`K` is the strike.
+
+    See :func:`pfhedge.nn.functional.realized_variance`
+    for the definition of the realized variance.
 
-    where :math:`\sigma^2` is the realized variance of the underlying asset
-    until maturity and :math:`K` is the strike variance (``strike``).
-    See :func:`pfhedge.nn.functional.realized_variance` for the definition of
-    the realized variance.
+    .. seealso::
+        :func:`pfhedge.nn.functional.realized_variance`
 
     Args:
         underlier (:class:`BasePrimary`): The underlying instrument.
diff --git a/pfhedge/instruments/primary/base.py b/pfhedge/instruments/primary/base.py
index 602d6b93..aafd05ad 100644
--- a/pfhedge/instruments/primary/base.py
+++ b/pfhedge/instruments/primary/base.py
@@ -1,5 +1,6 @@
 from abc import abstractmethod
 from collections import OrderedDict
+from typing import Any
 from typing import Dict
 from typing import Iterator
 from typing import Optional
@@ -40,18 +41,18 @@ class BasePrimary(BaseInstrument):
 
     dt: float
     cost: float
-    _buffers: Dict[str, Optional[Tensor]]
+    _buffers: Dict[str, Tensor]
     dtype: Optional[torch.dtype]
     device: Optional[torch.device]
 
     def __init__(self) -> None:
         super().__init__()
         self._buffers = OrderedDict()
-        self.register_buffer("spot", None)
 
     @property
     def default_init_state(self) -> Tuple[TensorOrScalar, ...]:
         """Returns the default initial state of simulation."""
+        return ()
 
     # TODO(simaki): Remove @no_type_check once BrownianStock and HestonStock get
     #   unified signatures.
@@ -78,7 +79,7 @@ def simulate(
                 (See :attr:`default_init_state`).
         """
 
-    def register_buffer(self, name: str, tensor: Optional[Tensor]) -> None:
+    def register_buffer(self, name: str, tensor: Tensor) -> None:
         """Adds a buffer to the instrument.
         The dtype and device of the buffer are the instrument's dtype and device.
 
@@ -133,14 +134,22 @@ def buffers(self) -> Iterator[Tensor]:
         for _, buffer in self.named_buffers():
             yield buffer
 
-    def __getattr__(self, name: str) -> Tensor:
+    def get_buffer(self, name: str) -> Tensor:
+        """Returns the buffer given by target if it exists, otherwise throws an error.
+
+        Args:
+            name (str): the name of the buffer.
+
+        Returns:
+            torch.Tensor
+        """
         if "_buffers" in self.__dict__:
-            _buffers = self.__dict__["_buffers"]
-            if name in _buffers:
-                return _buffers[name]
-        raise AttributeError(
-            "'{}' object has no attribute '{}'".format(type(self).__name__, name)
-        )
+            if name in self._buffers:
+                return self._buffers[name]
+        raise AttributeError(self._get_name() + " has no buffer named " + name)
+
+    def __getattr__(self, name: str) -> Tensor:
+        return self.get_buffer(name)
 
     @property
     def spot(self) -> Tensor:
@@ -158,7 +167,7 @@ def spot(self) -> Tensor:
     def is_listed(self) -> bool:
         return True
 
-    def to(self: T, *args, **kwargs) -> T:
+    def to(self: T, *args: Any, **kwargs: Any) -> T:
         device, dtype, *_ = self._parse_to(*args, **kwargs)
 
         if dtype is not None and not dtype.is_floating_point:
@@ -178,7 +187,7 @@ def to(self: T, *args, **kwargs) -> T:
         return self
 
     @staticmethod
-    def _parse_to(*args, **kwargs):
+    def _parse_to(*args: Any, **kwargs: Any):
         # Can be called as:
         #   to(device=None, dtype=None)
         #   to(tensor)
diff --git a/pfhedge/nn/functional.py b/pfhedge/nn/functional.py
index 6b2becb9..8b81fd83 100644
--- a/pfhedge/nn/functional.py
+++ b/pfhedge/nn/functional.py
@@ -8,12 +8,16 @@
 from torch.distributions.normal import Normal
 from torch.distributions.utils import broadcast_all
 
+import pfhedge.autogreek as autogreek
 from pfhedge._utils.typing import TensorOrScalar
 
 
 def european_payoff(input: Tensor, call: bool = True, strike: float = 1.0) -> Tensor:
     """Returns the payoff of a European option.
 
+    .. seealso::
+        - :class:`pfhedge.instruments.EuropeanOption`
+
     Args:
         input (torch.Tensor): The input tensor representing the price trajectory.
         call (bool, default=True): Specifies whether the option is call or put.
@@ -37,6 +41,9 @@ def european_payoff(input: Tensor, call: bool = True, strike: float = 1.0) -> Te
 def lookback_payoff(input: Tensor, call: bool = True, strike: float = 1.0) -> Tensor:
     """Returns the payoff of a lookback option with a fixed strike.
 
+    .. seealso::
+        - :class:`pfhedge.instruments.LookbackOption`
+
     Args:
         input (torch.Tensor): The input tensor representing the price trajectory.
         call (bool, default=True): Specifies whether the option is call or put.
@@ -62,6 +69,9 @@ def american_binary_payoff(
 ) -> Tensor:
     """Returns the payoff of an American binary option.
 
+    .. seealso::
+        - :class:`pfhedge.instruments.AmericanBinaryOption`
+
     Args:
         input (torch.Tensor): The input tensor representing the price trajectory.
         call (bool, default=True): Specifies whether the option is call or put.
@@ -87,6 +97,9 @@ def european_binary_payoff(
 ) -> Tensor:
     """Returns the payoff of a European binary option.
 
+    .. seealso::
+        - :class:`pfhedge.instruments.EuropeanBinaryOption`
+
     Args:
         input (torch.Tensor): The input tensor representing the price trajectory.
         call (bool, default=True): Specifies whether the option is call or put.
@@ -108,14 +121,18 @@ def european_binary_payoff(
 
 
 def european_forward_start_payoff(
-    input: Tensor, strike: float = 1.0, start_index: int = 0
+    input: Tensor, strike: float = 1.0, start_index: int = 0, end_index: int = -1
 ) -> Tensor:
     """Returns the payoff of a European forward start option.
 
+    .. seealso::
+        - :class:`pfhedge.instruments.EuropeanForwardStartOption`
+
     Args:
         input (torch.Tensor): The input tensor representing the price trajectory.
-        start_index (torch.Tensor): The time index at which the option starts.
         strike (float, default=1.0): The strike price of the option.
+        start_index (int, default=0): The time index at which the option starts.
+        end_index (int, default=-1): The time index at which the option ends.
 
     Shape:
         - input: :math:`(*, T)` where
@@ -126,7 +143,7 @@ def european_forward_start_payoff(
     Returns:
         torch.Tensor
     """
-    return fn.relu(input[..., -1] / input[..., start_index] - strike)
+    return fn.relu(input[..., end_index] / input[..., start_index] - strike)
 
 
 def exp_utility(input: Tensor, a: float = 1.0) -> Tensor:
@@ -517,11 +534,15 @@ def npdf(input: Tensor) -> Tensor:
     return Normal(0.0, 1.0).log_prob(input).exp()
 
 
-def d1(log_moneyness: Tensor, time_to_maturity: Tensor, volatility: Tensor) -> Tensor:
+def d1(
+    log_moneyness: TensorOrScalar,
+    time_to_maturity: TensorOrScalar,
+    volatility: TensorOrScalar,
+) -> Tensor:
     r"""Returns :math:`d_1` in the Black-Scholes formula.
 
     .. math::
-        d_1 = \frac{s + \frac12 \sigma^2 t}{\sigma \sqrt{t}}
+        d_1 = \frac{s}{\sigma \sqrt{t}} + \frac{\sigma \sqrt{t}}{2}
 
     where
     :math:`s` is the log moneyness,
@@ -544,19 +565,21 @@ def d1(log_moneyness: Tensor, time_to_maturity: Tensor, volatility: Tensor) -> T
         raise ValueError("all elements in time_to_maturity have to be non-negative")
     if not (v >= 0).all():
         raise ValueError("all elements in volatility have to be non-negative")
-    numerator = s + (v.square() / 2) * t
-    denominator = v * t.sqrt()
-    output = numerator / denominator
-    return torch.where(
-        (numerator == 0).logical_and(denominator == 0), torch.zeros_like(output), output
-    )
+    variance = v * t.sqrt()
+    output = s / variance + variance / 2
+    # TODO(simaki): Replace zeros_like with 0.0 once https://github.com/pytorch/pytorch/pull/62084 is merged
+    return output.where((s != 0).logical_or(variance != 0), torch.zeros_like(output))
 
 
-def d2(log_moneyness: Tensor, time_to_maturity: Tensor, volatility: Tensor) -> Tensor:
+def d2(
+    log_moneyness: TensorOrScalar,
+    time_to_maturity: TensorOrScalar,
+    volatility: TensorOrScalar,
+) -> Tensor:
     r"""Returns :math:`d_2` in the Black-Scholes formula.
 
     .. math::
-        d_2 = \frac{s - \frac12 \sigma^2 t}{\sigma \sqrt{t}}
+        d_2 = \frac{s}{\sigma \sqrt{t}} - \frac{\sigma \sqrt{t}}{2}
 
     where
     :math:`s` is the log moneyness,
@@ -579,12 +602,10 @@ def d2(log_moneyness: Tensor, time_to_maturity: Tensor, volatility: Tensor) -> T
         raise ValueError("all elements in time_to_maturity have to be non-negative")
     if not (v >= 0).all():
         raise ValueError("all elements in volatility have to be non-negative")
-    numerator = s - (v.square() / 2) * t
-    denominator = v * t.sqrt()
-    output = numerator / denominator
-    return torch.where(
-        (numerator == 0).logical_and(denominator == 0), torch.zeros_like(output), output
-    )
+    variance = v * t.sqrt()
+    output = s / variance - variance / 2
+    # TODO(simaki): Replace zeros_like with 0.0 once https://github.com/pytorch/pytorch/pull/62084 is merged
+    return output.where((s != 0).logical_or(variance != 0), torch.zeros_like(output))
 
 
 def ww_width(
@@ -633,3 +654,474 @@ def svi_variance(
     """
     k_m = torch.as_tensor(input - m)  # k - m
     return a + b * (rho * k_m + (k_m.square() + sigma ** 2).sqrt())
+
+
+def bilerp(
+    input1: Tensor,
+    input2: Tensor,
+    input3: Tensor,
+    input4: Tensor,
+    weight1: TensorOrScalar,
+    weight2: TensorOrScalar,
+) -> Tensor:
+    r"""Does a bilinear interpolation of four tensors based on a scalar or tensor weights and
+    returns the resulting tensor.
+
+    The output is given by
+
+    .. math::
+        \text{output}_i
+        & = (1 - w_1) (1 - w_2) \cdot \text{input1}_i
+        + w_1 (1 - w_2) \cdot \text{input2}_i \\
+        & \quad + (1 - w_1) w_2 \cdot \text{input3}_i
+        + w_1 w_2 \cdot \text{input4}_i ,
+
+    where :math:`w_1` and :math:`w_2` are the weights.
+
+    The shapes of inputs must be broadcastable.
+    If ``weight`` is a tensor, then the shapes of ``weight`` must also be broadcastable.
+
+    Args:
+        input1 (torch.Tensor): The input tensor.
+        input2 (torch.Tensor): The input tensor.
+        input3 (torch.Tensor): The input tensor.
+        input4 (torch.Tensor): The input tensor.
+        weight1 (float or torch.Tensor): The weight tensor.
+        weight2 (float or torch.Tensor): The weight tensor.
+
+    Returns:
+        torch.Tensor
+    """
+    lerp1 = torch.lerp(input1, input2, weight1)
+    lerp2 = torch.lerp(input3, input4, weight1)
+    return torch.lerp(lerp1, lerp2, weight2)
+
+
+def bs_european_price(
+    log_moneyness: Tensor,
+    time_to_maturity: Tensor,
+    volatility: Tensor,
+    strike: TensorOrScalar = 1.0,
+    call: bool = True,
+) -> Tensor:
+    """Returns Black-Scholes price of a European option.
+
+    See :func:`pfhedge.nn.BSEuropeanOption.price` for details.
+    """
+    s, t, v = broadcast_all(log_moneyness, time_to_maturity, volatility)
+
+    price = strike * (s.exp() * ncdf(d1(s, t, v)) - ncdf(d2(s, t, v)))
+    price = price + strike * (1 - s.exp()) if not call else price  # put-call parity
+
+    return price
+
+
+def bs_european_delta(
+    log_moneyness: Tensor,
+    time_to_maturity: Tensor,
+    volatility: Tensor,
+    call: bool = True,
+) -> Tensor:
+    """Returns Black-Scholes delta of a European option.
+
+    See :func:`pfhedge.nn.BSEuropeanOption.delta` for details.
+    """
+    s, t, v = broadcast_all(log_moneyness, time_to_maturity, volatility)
+
+    delta = ncdf(d1(s, t, v))
+    delta = delta - 1 if not call else delta  # put-call parity
+
+    return delta
+
+
+def bs_european_gamma(
+    log_moneyness: Tensor,
+    time_to_maturity: Tensor,
+    volatility: Tensor,
+    strike: TensorOrScalar = 1.0,
+) -> Tensor:
+    """Returns Black-Scholes gamma of a European option.
+
+    See :func:`pfhedge.nn.BSEuropeanOption.gamma` for details.
+    """
+    s, t, v = broadcast_all(log_moneyness, time_to_maturity, volatility)
+    spot = strike * s.exp()
+    numerator = npdf(d1(s, t, v))
+    denominator = spot * v * t.sqrt()
+    output = numerator / denominator
+    return torch.where(
+        (numerator == 0).logical_and(denominator == 0), torch.zeros_like(output), output
+    )
+
+
+def bs_european_vega(
+    log_moneyness: Tensor,
+    time_to_maturity: Tensor,
+    volatility: Tensor,
+    strike: TensorOrScalar,
+) -> Tensor:
+    """Returns Black-Scholes vega of a European option.
+
+    See :func:`pfhedge.nn.BSEuropeanOption.vega` for details.
+    """
+    s, t, v = broadcast_all(log_moneyness, time_to_maturity, volatility)
+    price = strike * s.exp()
+    return npdf(d1(s, t, v)) * price * t.sqrt()
+
+
+def bs_european_theta(
+    log_moneyness: Tensor,
+    time_to_maturity: Tensor,
+    volatility: Tensor,
+    strike: TensorOrScalar,
+) -> Tensor:
+    """Returns Black-Scholes theta of a European option.
+
+    See :func:`pfhedge.nn.BSEuropeanOption.theta` for details.
+    """
+    s, t, v = broadcast_all(log_moneyness, time_to_maturity, volatility)
+    price = strike * s.exp()
+    numerator = -npdf(d1(s, t, v)) * price * v
+    denominator = 2 * t.sqrt()
+    output = numerator / denominator
+    return torch.where(
+        (numerator == 0).logical_and(denominator == 0), torch.zeros_like(output), output
+    )
+
+
+def bs_european_binary_price(
+    log_moneyness: Tensor,
+    time_to_maturity: Tensor,
+    volatility: Tensor,
+    call: bool = True,
+) -> Tensor:
+    """Returns Black-Scholes price of a European binary option.
+
+    See :func:`pfhedge.nn.BSEuropeanBinaryOption.price` for details.
+    """
+    s, t, v = broadcast_all(log_moneyness, time_to_maturity, volatility)
+
+    price = ncdf(d2(s, t, v))
+    price = 1.0 - price if not call else price  # put-call parity
+
+    return price
+
+
+def bs_european_binary_delta(
+    log_moneyness: Tensor,
+    time_to_maturity: Tensor,
+    volatility: Tensor,
+    call: bool = True,
+    strike: TensorOrScalar = 1.0,
+) -> Tensor:
+    """Returns Black-Scholes delta of a European binary option.
+
+    See :func:`pfhedge.nn.BSEuropeanBinaryOption.delta` for details.
+    """
+    s, t, v = broadcast_all(log_moneyness, time_to_maturity, volatility)
+
+    spot = s.exp() * strike
+
+    numerator = npdf(d2(s, t, v))
+    denominator = spot * v * t.sqrt()
+    delta = numerator / denominator
+    delta = torch.where(
+        (numerator == 0).logical_and(denominator == 0), torch.zeros_like(delta), delta
+    )
+    delta = -delta if not call else delta  # put-call parity
+
+    return delta
+
+
+def bs_european_binary_gamma(
+    log_moneyness: Tensor,
+    time_to_maturity: Tensor,
+    volatility: Tensor,
+    strike: TensorOrScalar = 1.0,
+) -> Tensor:
+    """Returns Black-Scholes gamma of a European binary option.
+
+    See :func:`pfhedge.nn.BSEuropeanBinaryOption.gamma` for details.
+    """
+    # TODO(simaki): Directly compute gamma.
+    return autogreek.gamma(
+        bs_european_binary_price,
+        log_moneyness=log_moneyness,
+        time_to_maturity=time_to_maturity,
+        volatility=volatility,
+        strike=strike,
+    )
+
+
+def bs_european_binary_vega(
+    log_moneyness: Tensor,
+    time_to_maturity: Tensor,
+    volatility: Tensor,
+    strike: TensorOrScalar = 1.0,
+) -> Tensor:
+    """Returns Black-Scholes vega of a European binary option.
+
+    See :func:`pfhedge.nn.BSEuropeanBinaryOption.vega` for details.
+    """
+    # TODO(simaki): Directly compute gamma.
+    return autogreek.vega(
+        bs_european_binary_price,
+        log_moneyness=log_moneyness,
+        time_to_maturity=time_to_maturity,
+        volatility=volatility,
+        strike=strike,
+    )
+
+
+def bs_european_binary_theta(
+    log_moneyness: Tensor,
+    time_to_maturity: Tensor,
+    volatility: Tensor,
+    strike: TensorOrScalar = 1.0,
+) -> Tensor:
+    """Returns Black-Scholes theta of a European binary option.
+
+    See :func:`pfhedge.nn.BSEuropeanBinaryOption.theta` for details.
+    """
+    # TODO(simaki): Directly compute theta.
+    return autogreek.theta(
+        bs_european_binary_price,
+        log_moneyness=log_moneyness,
+        time_to_maturity=time_to_maturity,
+        volatility=volatility,
+        strike=strike,
+    )
+
+
+def bs_american_binary_price(
+    log_moneyness: Tensor,
+    max_log_moneyness: Tensor,
+    time_to_maturity: Tensor,
+    volatility: Tensor,
+) -> Tensor:
+    """Returns Black-Scholes price of an American binary option.
+
+    See :func:`pfhedge.nn.BSAmericanBinaryOption.price` for details.
+    """
+    # This formula is derived using the results in Section 7.3.3 of Shreve's book.
+    # Price is I_2 - I_4 where the interval of integration is [k --> -inf, b].
+    # By this substitution we get N([log(S(0) / K) + ...] / sigma T) --> 1.
+
+    s, t, v = broadcast_all(log_moneyness, time_to_maturity, volatility)
+    p = ncdf(d2(s, t, v)) + s.exp() * (1 - ncdf(d2(-s, t, v)))
+
+    return p.where(max_log_moneyness < 0, torch.ones_like(p))
+
+
+def bs_american_binary_delta(
+    log_moneyness: Tensor,
+    max_log_moneyness: Tensor,
+    time_to_maturity: Tensor,
+    volatility: Tensor,
+    strike: TensorOrScalar,
+) -> Tensor:
+    """Returns Black-Scholes delta of an American binary option.
+
+    See :func:`pfhedge.nn.BSAmericanBinaryOption.delta` for details.
+    """
+    s, t, v = broadcast_all(log_moneyness, time_to_maturity, volatility)
+    spot = s.exp() * strike
+    # ToDo: fix 0/0 issue
+    p = (
+        npdf(d2(s, t, v)) / (spot * v * t.sqrt())
+        - (1 - ncdf(d2(-s, t, v))) / strike
+        + npdf(d2(-s, t, v)) / (strike * v * t.sqrt())
+    )
+
+    return p.where(max_log_moneyness < 0, torch.zeros_like(p))
+
+
+def bs_american_binary_gamma(
+    log_moneyness: Tensor,
+    max_log_moneyness: Tensor,
+    time_to_maturity: Tensor,
+    volatility: Tensor,
+    strike: TensorOrScalar,
+) -> Tensor:
+    """Returns Black-Scholes gamma of an American binary option.
+
+    See :func:`pfhedge.nn.BSAmericanBinaryOption.gamma` for details.
+    """
+    # TODO(simaki): Compute analytically
+    return autogreek.gamma(
+        bs_american_binary_price,
+        log_moneyness=log_moneyness,
+        max_log_moneyness=max_log_moneyness,
+        time_to_maturity=time_to_maturity,
+        volatility=volatility,
+        strike=strike,
+    )
+
+
+def bs_american_binary_vega(
+    log_moneyness: Tensor,
+    max_log_moneyness: Tensor,
+    time_to_maturity: Tensor,
+    volatility: Tensor,
+    strike: TensorOrScalar,
+) -> Tensor:
+    """Returns Black-Scholes vega of an American binary option.
+
+    See :func:`pfhedge.nn.BSAmericanBinaryOption.vega` for details.
+    """
+    # TODO(simaki): Compute analytically
+    return autogreek.vega(
+        bs_american_binary_price,
+        log_moneyness=log_moneyness,
+        max_log_moneyness=max_log_moneyness,
+        time_to_maturity=time_to_maturity,
+        volatility=volatility,
+        strike=strike,
+    )
+
+
+def bs_american_binary_theta(
+    log_moneyness: Tensor,
+    max_log_moneyness: Tensor,
+    time_to_maturity: Tensor,
+    volatility: Tensor,
+    strike: TensorOrScalar,
+) -> Tensor:
+    """Returns Black-Scholes theta of an American binary option.
+
+    See :func:`pfhedge.nn.BSAmericanBinaryOption.theta` for details.
+    """
+    # TODO(simaki): Compute analytically
+    return autogreek.theta(
+        bs_american_binary_price,
+        log_moneyness=log_moneyness,
+        max_log_moneyness=max_log_moneyness,
+        time_to_maturity=time_to_maturity,
+        volatility=volatility,
+        strike=strike,
+    )
+
+
+def bs_lookback_price(
+    log_moneyness: Tensor,
+    max_log_moneyness: Tensor,
+    time_to_maturity: Tensor,
+    volatility: Tensor,
+    strike: TensorOrScalar,
+) -> Tensor:
+    """Returns Black-Scholes price of a lookback option.
+
+    See :func:`pfhedge.nn.BSLookbackOption.price` for details.
+    """
+    s, m, t, v = map(
+        torch.as_tensor,
+        (log_moneyness, max_log_moneyness, time_to_maturity, volatility),
+    )
+
+    spot = s.exp() * strike
+    max = m.exp() * strike
+    d1_value = d1(s, t, v)
+    d2_value = d2(s, t, v)
+    m1 = d1(s - m, t, v)  # d' in the paper
+    m2 = d2(s - m, t, v)
+
+    # when max < strike
+    price_0 = spot * (
+        ncdf(d1_value) + v * t.sqrt() * (d1_value * ncdf(d1_value) + npdf(d1_value))
+    ) - strike * ncdf(d2_value)
+    # when max >= strike
+    price_1 = (
+        spot * (ncdf(m1) + v * t.sqrt() * (m1 * ncdf(m1) + npdf(m1)))
+        - strike
+        + max * (1 - ncdf(m2))
+    )
+
+    return torch.where(max < strike, price_0, price_1)
+
+
+def bs_lookback_delta(
+    log_moneyness: Tensor,
+    max_log_moneyness: Tensor,
+    time_to_maturity: Tensor,
+    volatility: Tensor,
+    strike: TensorOrScalar,
+) -> Tensor:
+    """Returns Black-Scholes delta of a lookback option.
+
+    See :func:`pfhedge.nn.BSLookbackOption.delta` for details.
+    """
+    # TODO(simaki): Calculate analytically
+    return autogreek.delta(
+        bs_lookback_price,
+        log_moneyness=log_moneyness,
+        max_log_moneyness=max_log_moneyness,
+        time_to_maturity=time_to_maturity,
+        volatility=volatility,
+        strike=strike,
+    )
+
+
+def bs_lookback_gamma(
+    log_moneyness: Tensor,
+    max_log_moneyness: Tensor,
+    time_to_maturity: Tensor,
+    volatility: Tensor,
+    strike: TensorOrScalar,
+) -> Tensor:
+    """Returns Black-Scholes gamma of a lookback option.
+
+    See :func:`pfhedge.nn.BSLookbackOption.gamma` for details.
+    """
+    # TODO(simaki): Calculate analytically
+    return autogreek.gamma(
+        bs_lookback_price,
+        log_moneyness=log_moneyness,
+        max_log_moneyness=max_log_moneyness,
+        time_to_maturity=time_to_maturity,
+        volatility=volatility,
+        strike=strike,
+    )
+
+
+def bs_lookback_vega(
+    log_moneyness: Tensor,
+    max_log_moneyness: Tensor,
+    time_to_maturity: Tensor,
+    volatility: Tensor,
+    strike: TensorOrScalar,
+) -> Tensor:
+    """Returns Black-Scholes vega of a lookback option.
+
+    See :func:`pfhedge.nn.BSLookbackOption.vega` for details.
+    """
+    # TODO(simaki): Calculate analytically
+    return autogreek.vega(
+        bs_lookback_price,
+        log_moneyness=log_moneyness,
+        max_log_moneyness=max_log_moneyness,
+        time_to_maturity=time_to_maturity,
+        volatility=volatility,
+        strike=strike,
+    )
+
+
+def bs_lookback_theta(
+    log_moneyness: Tensor,
+    max_log_moneyness: Tensor,
+    time_to_maturity: Tensor,
+    volatility: Tensor,
+    strike: TensorOrScalar,
+) -> Tensor:
+    """Returns Black-Scholes theta of a lookback option.
+
+    See :func:`pfhedge.nn.BSLookbackOption.theta` for details.
+    """
+    # TODO(simaki): Calculate analytically
+    return autogreek.theta(
+        bs_lookback_price,
+        log_moneyness=log_moneyness,
+        max_log_moneyness=max_log_moneyness,
+        time_to_maturity=time_to_maturity,
+        volatility=volatility,
+        strike=strike,
+    )
diff --git a/pfhedge/nn/modules/bs/_base.py b/pfhedge/nn/modules/bs/_base.py
index c9e15c3e..84406766 100644
--- a/pfhedge/nn/modules/bs/_base.py
+++ b/pfhedge/nn/modules/bs/_base.py
@@ -1,5 +1,8 @@
 from inspect import signature
 from typing import List
+from typing import Optional
+from typing import Tuple
+from typing import Union
 from typing import no_type_check
 
 import torch
@@ -7,6 +10,10 @@
 from torch.nn import Module
 
 import pfhedge.autogreek as autogreek
+from pfhedge.instruments import AmericanBinaryOption
+from pfhedge.instruments import EuropeanBinaryOption
+from pfhedge.instruments import EuropeanOption
+from pfhedge.instruments import LookbackOption
 
 
 class BSModuleMixin(Module):
@@ -85,3 +92,82 @@ def inputs(self) -> List[str]:
             list
         """
         return list(signature(self.delta).parameters.keys())
+
+
+def acquire_params_from_derivative_0(
+    derivative: Optional[
+        Union[
+            EuropeanOption, EuropeanBinaryOption, AmericanBinaryOption, LookbackOption
+        ]
+    ],
+    log_moneyness: Optional[Tensor] = None,
+    time_to_maturity: Optional[Tensor] = None,
+) -> Tuple[Tensor, Tensor]:
+    if log_moneyness is None:
+        if derivative is None:
+            raise ValueError(
+                "log_moneyness is required if derivative is not set at this initialization."
+            )
+        if derivative.ul().spot is None:
+            raise AttributeError("please simulate first")
+        log_moneyness = derivative.log_moneyness()
+    if time_to_maturity is None:
+        if derivative is None:
+            raise ValueError(
+                "time_to_maturity is required if derivative is not set at this initialization."
+            )
+        time_to_maturity = derivative.time_to_maturity()
+    return log_moneyness, time_to_maturity
+
+
+def acquire_params_from_derivative_1(
+    derivative: Optional[
+        Union[
+            EuropeanOption, EuropeanBinaryOption, AmericanBinaryOption, LookbackOption
+        ]
+    ],
+    log_moneyness: Optional[Tensor] = None,
+    time_to_maturity: Optional[Tensor] = None,
+    volatility: Optional[Tensor] = None,
+) -> Tuple[Tensor, Tensor, Tensor]:
+    log_moneyness, time_to_maturity = acquire_params_from_derivative_0(
+        derivative=derivative,
+        log_moneyness=log_moneyness,
+        time_to_maturity=time_to_maturity,
+    )
+    if volatility is None:
+        if derivative is None:
+            raise ValueError(
+                "time_to_maturity is required if derivative is not set at this initialization."
+            )
+        if derivative.ul().volatility is None:
+            raise AttributeError(
+                "please simulate first and check if volatility exists in the derivative's underlier."
+            )
+        volatility = derivative.ul().volatility
+    return log_moneyness, time_to_maturity, volatility
+
+
+def acquire_params_from_derivative_2(
+    derivative: Optional[Union[AmericanBinaryOption, LookbackOption]],
+    log_moneyness: Optional[Tensor] = None,
+    max_log_moneyness: Optional[Tensor] = None,
+    time_to_maturity: Optional[Tensor] = None,
+    volatility: Optional[Tensor] = None,
+) -> Tuple[Tensor, Tensor, Tensor, Tensor]:
+    log_moneyness, time_to_maturity, volatility = acquire_params_from_derivative_1(
+        derivative=derivative,
+        log_moneyness=log_moneyness,
+        time_to_maturity=time_to_maturity,
+        volatility=volatility,
+    )
+    if max_log_moneyness is None:
+        if derivative is None:
+            raise ValueError(
+                "max_log_moneyness is required if derivative is not set at this initialization."
+            )
+        if derivative.ul().spot is None:
+            raise AttributeError("please simulate first")
+        max_log_moneyness = derivative.max_log_moneyness()
+
+    return log_moneyness, max_log_moneyness, time_to_maturity, volatility
diff --git a/pfhedge/nn/modules/bs/american_binary.py b/pfhedge/nn/modules/bs/american_binary.py
index 5b06a1f6..86b73647 100644
--- a/pfhedge/nn/modules/bs/american_binary.py
+++ b/pfhedge/nn/modules/bs/american_binary.py
@@ -1,25 +1,30 @@
 from math import sqrt
+from typing import Optional
 
 import torch
 from torch import Tensor
-from torch.distributions.utils import broadcast_all
 
 from pfhedge._utils.bisect import find_implied_volatility
 from pfhedge._utils.doc import _set_attr_and_docstring
 from pfhedge._utils.str import _format_float
-from pfhedge.nn.functional import d1
-from pfhedge.nn.functional import d2
-from pfhedge.nn.functional import ncdf
-from pfhedge.nn.functional import npdf
+from pfhedge.instruments import AmericanBinaryOption
+from pfhedge.nn.functional import bs_american_binary_delta
+from pfhedge.nn.functional import bs_american_binary_price
 
 from ._base import BSModuleMixin
+from ._base import acquire_params_from_derivative_0
+from ._base import acquire_params_from_derivative_2
+from .black_scholes import BlackScholesModuleFactory
 
 
 class BSAmericanBinaryOption(BSModuleMixin):
-    """Black-Scholes formula for an american binary option.
+    """Black-Scholes formula for an American binary option.
 
     Note:
-        Risk-free rate is set to zero.
+        - The formulas are for continuous monitoring while
+          :class:`pfhedge.instruments.AmericanBinaryOption` monitors spot prices discretely.
+          To get adjustment for discrete monitoring, see, for instance,
+          Broadie, Glasserman, and Kou (1999).
 
     Args:
         call (bool, default=True): Specifies whether the option is call or put.
@@ -41,6 +46,9 @@ class BSAmericanBinaryOption(BSModuleMixin):
     References:
         - Shreve, S.E., 2004. Stochastic calculus for finance II:
           Continuous-time models (Vol. 11). Springer Science & Business Media.
+        - Broadie, M., Glasserman, P. and Kou, S.G., 1999.
+          Connecting discrete and continuous path-dependent options.
+          Finance and Stochastics, 3(1), pp.55-82.
 
     Examples:
         >>> from pfhedge.nn import BSAmericanBinaryOption
@@ -58,7 +66,12 @@ class BSAmericanBinaryOption(BSModuleMixin):
                 [...]])
     """
 
-    def __init__(self, call: bool = True, strike: float = 1.0):
+    def __init__(
+        self,
+        call: bool = True,
+        strike: float = 1.0,
+        derivative: Optional[AmericanBinaryOption] = None,
+    ):
         if not call:
             raise ValueError(
                 f"{self.__class__.__name__} for a put option is not yet supported."
@@ -67,6 +80,7 @@ def __init__(self, call: bool = True, strike: float = 1.0):
         super().__init__()
         self.call = call
         self.strike = strike
+        self.derivative = derivative
 
     @classmethod
     def from_derivative(cls, derivative):
@@ -88,7 +102,9 @@ def from_derivative(cls, derivative):
             >>> m
             BSAmericanBinaryOption(strike=1.1000)
         """
-        return cls(call=derivative.call, strike=derivative.strike)
+        return cls(
+            call=derivative.call, strike=derivative.strike, derivative=derivative
+        )
 
     def extra_repr(self) -> str:
         params = []
@@ -97,18 +113,18 @@ def extra_repr(self) -> str:
 
     def price(
         self,
-        log_moneyness: Tensor,
-        max_log_moneyness: Tensor,
-        time_to_maturity: Tensor,
-        volatility: Tensor,
+        log_moneyness: Optional[Tensor] = None,
+        max_log_moneyness: Optional[Tensor] = None,
+        time_to_maturity: Optional[Tensor] = None,
+        volatility: Optional[Tensor] = None,
     ) -> Tensor:
         """Returns price of the derivative.
 
         Args:
-            log_moneyness (torch.Tensor): Log moneyness of the underlying asset.
-            max_log_moneyness (torch.Tensor): Cumulative maximum of the log moneyness.
-            time_to_maturity (torch.Tensor): Time to expiry of the option.
-            volatility (torch.Tensor): Volatility of the underlying asset.
+            log_moneyness (torch.Tensor, optional): Log moneyness of the underlying asset.
+            max_log_moneyness (torch.Tensor, optional): Cumulative maximum of the log moneyness.
+            time_to_maturity (torch.Tensor, optional): Time to expiry of the option.
+            volatility (torch.Tensor, optional): Volatility of the underlying asset.
 
         Shape:
             - log_moneyness: :math:`(N, *)` where
@@ -120,30 +136,46 @@ def price(
 
         Returns:
             torch.Tensor
+
+        Note:
+            Parameters are not optional if the module has not accepted a derivative in its initialization.
         """
         # This formula is derived using the results in Section 7.3.3 of Shreve's book.
         # Price is I_2 - I_4 where the interval of integration is [k --> -inf, b].
         # By this substitution we get N([log(S(0) / K) + ...] / sigma T) --> 1.
-
-        s, t, v = broadcast_all(log_moneyness, time_to_maturity, volatility)
-        p = ncdf(d2(s, t, v)) + s.exp() * (1 - ncdf(d2(-s, t, v)))
-
-        return p.where(max_log_moneyness < 0, torch.ones_like(p))
+        (
+            log_moneyness,
+            max_log_moneyness,
+            time_to_maturity,
+            volatility,
+        ) = acquire_params_from_derivative_2(
+            self.derivative,
+            log_moneyness,
+            max_log_moneyness,
+            time_to_maturity,
+            volatility,
+        )
+        return bs_american_binary_price(
+            log_moneyness=log_moneyness,
+            max_log_moneyness=max_log_moneyness,
+            time_to_maturity=time_to_maturity,
+            volatility=volatility,
+        )
 
     def delta(
         self,
-        log_moneyness: Tensor,
-        max_log_moneyness: Tensor,
-        time_to_maturity: Tensor,
-        volatility: Tensor,
+        log_moneyness: Optional[Tensor] = None,
+        max_log_moneyness: Optional[Tensor] = None,
+        time_to_maturity: Optional[Tensor] = None,
+        volatility: Optional[Tensor] = None,
     ) -> Tensor:
         """Returns delta of the derivative.
 
         Args:
-            log_moneyness (torch.Tensor): Log moneyness of the underlying asset.
-            max_log_moneyness (torch.Tensor): Cumulative maximum of the log moneyness.
-            time_to_maturity (torch.Tensor): Time to expiry of the option.
-            volatility (torch.Tensor): Volatility of the underlying asset.
+            log_moneyness (torch.Tensor, optional): Log moneyness of the underlying asset.
+            max_log_moneyness (torch.Tensor, optional): Cumulative maximum of the log moneyness.
+            time_to_maturity (torch.Tensor, optional): Time to expiry of the option.
+            volatility (torch.Tensor, optional): Volatility of the underlying asset.
 
         Shape:
             - log_moneyness: :math:`(N, *)` where
@@ -155,33 +187,45 @@ def delta(
 
         Returns:
             torch.Tensor
+
+        Note:
+            Parameters are not optional if the module has not accepted a derivative in its initialization.
         """
-        s, t, v = broadcast_all(log_moneyness, time_to_maturity, volatility)
-        spor = s.exp() * self.strike
-        # ToDo: fix 0/0 issue
-        p = (
-            npdf(d2(s, t, v)) / (spor * v * t.sqrt())
-            - (1 - ncdf(d2(-s, t, v))) / self.strike
-            + npdf(d2(-s, t, v)) / (self.strike * v * t.sqrt())
+        (
+            log_moneyness,
+            max_log_moneyness,
+            time_to_maturity,
+            volatility,
+        ) = acquire_params_from_derivative_2(
+            self.derivative,
+            log_moneyness,
+            max_log_moneyness,
+            time_to_maturity,
+            volatility,
+        )
+        return bs_american_binary_delta(
+            log_moneyness=log_moneyness,
+            max_log_moneyness=max_log_moneyness,
+            time_to_maturity=time_to_maturity,
+            volatility=volatility,
+            strike=self.strike,
         )
-
-        return p.where(max_log_moneyness < 0, torch.zeros_like(p))
 
     @torch.enable_grad()
     def gamma(
         self,
-        log_moneyness: Tensor,
-        max_log_moneyness: Tensor,
-        time_to_maturity: Tensor,
-        volatility: Tensor,
+        log_moneyness: Optional[Tensor] = None,
+        max_log_moneyness: Optional[Tensor] = None,
+        time_to_maturity: Optional[Tensor] = None,
+        volatility: Optional[Tensor] = None,
     ) -> Tensor:
         """Returns gamma of the derivative.
 
         Args:
-            log_moneyness (torch.Tensor): Log moneyness of the underlying asset.
-            max_log_moneyness (torch.Tensor): Cumulative maximum of the log moneyness.
-            time_to_maturity (torch.Tensor): Time to expiry of the option.
-            volatility (torch.Tensor): Volatility of the underlying asset.
+            log_moneyness (torch.Tensor, optional): Log moneyness of the underlying asset.
+            max_log_moneyness (torch.Tensor, optional): Cumulative maximum of the log moneyness.
+            time_to_maturity (torch.Tensor, optional): Time to expiry of the option.
+            volatility (torch.Tensor, optional): Volatility of the underlying asset.
 
         Shape:
             - log_moneyness: :math:`(N, *)` where
@@ -193,30 +237,45 @@ def gamma(
 
         Returns:
             torch.Tensor
+
+        Note:
+            args are not optional if it doesn't accept derivative in this initialization.
         """
+        (
+            log_moneyness,
+            max_log_moneyness,
+            time_to_maturity,
+            volatility,
+        ) = acquire_params_from_derivative_2(
+            self.derivative,
+            log_moneyness,
+            max_log_moneyness,
+            time_to_maturity,
+            volatility,
+        )
         return super().gamma(
-            strike=self.strike,
             log_moneyness=log_moneyness,
             max_log_moneyness=max_log_moneyness,
             time_to_maturity=time_to_maturity,
             volatility=volatility,
+            strike=self.strike,
         )
 
     @torch.enable_grad()
     def vega(
         self,
-        log_moneyness: Tensor,
-        max_log_moneyness: Tensor,
-        time_to_maturity: Tensor,
-        volatility: Tensor,
+        log_moneyness: Optional[Tensor] = None,
+        max_log_moneyness: Optional[Tensor] = None,
+        time_to_maturity: Optional[Tensor] = None,
+        volatility: Optional[Tensor] = None,
     ) -> Tensor:
         """Returns vega of the derivative.
 
         Args:
-            log_moneyness (torch.Tensor): Log moneyness of the underlying asset.
-            max_log_moneyness (torch.Tensor): Cumulative maximum of the log moneyness.
-            time_to_maturity (torch.Tensor): Time to expiry of the option.
-            volatility (torch.Tensor): Volatility of the underlying asset.
+            log_moneyness (torch.Tensor, optional): Log moneyness of the underlying asset.
+            max_log_moneyness (torch.Tensor, optional): Cumulative maximum of the log moneyness.
+            time_to_maturity (torch.Tensor, optional): Time to expiry of the option.
+            volatility (torch.Tensor, optional): Volatility of the underlying asset.
 
         Shape:
             - log_moneyness: :math:`(N, *)` where
@@ -228,30 +287,45 @@ def vega(
 
         Returns:
             torch.Tensor
+
+        Note:
+            args are not optional if it doesn't accept derivative in this initialization.
         """
+        (
+            log_moneyness,
+            max_log_moneyness,
+            time_to_maturity,
+            volatility,
+        ) = acquire_params_from_derivative_2(
+            self.derivative,
+            log_moneyness,
+            max_log_moneyness,
+            time_to_maturity,
+            volatility,
+        )
         return super().vega(
-            strike=self.strike,
             log_moneyness=log_moneyness,
             max_log_moneyness=max_log_moneyness,
             time_to_maturity=time_to_maturity,
             volatility=volatility,
+            strike=self.strike,
         )
 
     @torch.enable_grad()
     def theta(
         self,
-        log_moneyness: Tensor,
-        max_log_moneyness: Tensor,
-        time_to_maturity: Tensor,
-        volatility: Tensor,
+        log_moneyness: Optional[Tensor] = None,
+        max_log_moneyness: Optional[Tensor] = None,
+        time_to_maturity: Optional[Tensor] = None,
+        volatility: Optional[Tensor] = None,
     ) -> Tensor:
         """Returns theta of the derivative.
 
         Args:
-            log_moneyness (torch.Tensor): Log moneyness of the underlying asset.
-            max_log_moneyness (torch.Tensor): Cumulative maximum of the log moneyness.
-            time_to_maturity (torch.Tensor): Time to expiry of the option.
-            volatility (torch.Tensor): Volatility of the underlying asset.
+            log_moneyness (torch.Tensor, optional): Log moneyness of the underlying asset.
+            max_log_moneyness (torch.Tensor, optional): Cumulative maximum of the log moneyness.
+            time_to_maturity (torch.Tensor, optional): Time to expiry of the option.
+            volatility (torch.Tensor, optional): Volatility of the underlying asset.
 
         Shape:
             - log_moneyness: :math:`(N, *)` where
@@ -266,32 +340,47 @@ def theta(
 
         Returns:
             torch.Tensor
+
+        Note:
+            args are not optional if it doesn't accept derivative in this initialization.
         """
+        (
+            log_moneyness,
+            max_log_moneyness,
+            time_to_maturity,
+            volatility,
+        ) = acquire_params_from_derivative_2(
+            self.derivative,
+            log_moneyness,
+            max_log_moneyness,
+            time_to_maturity,
+            volatility,
+        )
         return super().theta(
-            strike=self.strike,
             log_moneyness=log_moneyness,
             max_log_moneyness=max_log_moneyness,
             time_to_maturity=time_to_maturity,
             volatility=volatility,
+            strike=self.strike,
         )
 
     def implied_volatility(
         self,
-        log_moneyness: Tensor,
-        max_log_moneyness: Tensor,
-        time_to_maturity: Tensor,
-        price: Tensor,
+        log_moneyness: Optional[Tensor] = None,
+        max_log_moneyness: Optional[Tensor] = None,
+        time_to_maturity: Optional[Tensor] = None,
+        price: Optional[Tensor] = None,
         precision: float = 1e-6,
     ) -> Tensor:
         """Returns implied volatility of the derivative.
 
         Args:
-            log_moneyness (torch.Tensor): Log moneyness of the underlying asset.
-            max_log_moneyness (torch.Tensor): Cumulative maximum of the log moneyness.
-            time_to_maturity (torch.Tensor): Time to expiry of the option.
-            volatility (torch.Tensor): Volatility of the underlying asset.
-            precision (float, default=1e-6): Computational precision of the implied
-                volatility.
+            log_moneyness (torch.Tensor, optional): Log moneyness of the underlying asset.
+            max_log_moneyness (torch.Tensor, optional): Cumulative maximum of the log moneyness.
+            time_to_maturity (torch.Tensor, optional): Time to expiry of the option.
+            volatility (torch.Tensor, optional): Volatility of the underlying asset.
+            price (torch.Tensor): Price of the derivative.
+            precision (float, default=1e-6): Computational precision of the implied volatility.
 
         Shape:
             - log_moneyness: :math:`(N, *)` where
@@ -302,7 +391,18 @@ def implied_volatility(
 
         Returns:
             torch.Tensor
+
+        Note:
+            args are not optional if it doesn't accept derivative in this initialization.
+            price seems optional in typing, but it isn't. It is set for the compatibility to the previous versions.
         """
+        (log_moneyness, time_to_maturity) = acquire_params_from_derivative_0(
+            self.derivative, log_moneyness, time_to_maturity
+        )
+        if price is None:
+            raise ValueError(
+                "price is required in this method. None is set only for compatibility to the previous versions."
+            )
         return find_implied_volatility(
             self.price,
             price=price,
@@ -313,6 +413,9 @@ def implied_volatility(
         )
 
 
+factory = BlackScholesModuleFactory()
+factory.register_module("AmericanBinaryOption", BSAmericanBinaryOption)
+
 # Assign docstrings so they appear in Sphinx documentation
 _set_attr_and_docstring(BSAmericanBinaryOption, "inputs", BSModuleMixin.inputs)
 _set_attr_and_docstring(BSAmericanBinaryOption, "forward", BSModuleMixin.forward)
diff --git a/pfhedge/nn/modules/bs/black_scholes.py b/pfhedge/nn/modules/bs/black_scholes.py
index 27e68749..392b9c1f 100644
--- a/pfhedge/nn/modules/bs/black_scholes.py
+++ b/pfhedge/nn/modules/bs/black_scholes.py
@@ -1,13 +1,47 @@
+from collections import OrderedDict
 from typing import Callable
+from typing import Dict
+from typing import Iterator
 from typing import List
+from typing import Tuple
+from typing import Type
 
 from torch import Tensor
 from torch.nn import Module
 
-from .american_binary import BSAmericanBinaryOption
-from .european import BSEuropeanOption
-from .european_binary import BSEuropeanBinaryOption
-from .lookback import BSLookbackOption
+
+class BlackScholesModuleFactory:
+
+    _modules: Dict[str, Type[Module]]
+
+    # singleton
+    def __new__(cls, *args, **kwargs):
+        if not hasattr(cls, "_instance"):
+            cls._instance = super().__new__(cls)
+            cls._instance._modules = OrderedDict()
+        return cls._instance
+
+    def register_module(self, name: str, cls: Type[Module]) -> None:
+        self._modules[name] = cls
+
+    def named_modules(self) -> Iterator[Tuple[str, Type[Module]]]:
+        for name, module in self._modules.items():
+            if module is not None:
+                yield name, module
+
+    def names(self) -> Iterator[str]:
+        for name, _ in self.named_modules():
+            yield name
+
+    def features(self) -> Iterator[Type[Module]]:
+        for _, module in self.named_modules():
+            yield module
+
+    def get_class(self, name: str) -> Type[Module]:
+        return self._modules[name]
+
+    def get_class_from_derivative(self, derivative) -> Type[Module]:
+        return self.get_class(derivative.__class__.__name__).from_derivative(derivative)  # type: ignore
 
 
 class BlackScholes(Module):
@@ -42,16 +76,6 @@ class BlackScholes(Module):
         >>> m = BlackScholes(derivative)
         >>> m
         BSEuropeanOption(strike=1.1000)
-        >>> m.inputs()
-        ['log_moneyness', 'time_to_maturity', 'volatility']
-        >>> input = torch.tensor([
-        ...     [-0.01, 0.1, 0.2],
-        ...     [ 0.00, 0.1, 0.2],
-        ...     [ 0.01, 0.1, 0.2]])
-        >>> m(input)
-        tensor([[0.4497],
-                [0.5126],
-                [0.5752]])
 
         Instantiating :class:`BSLookbackOption` using a
         :class:`pfhedge.instruments.LookbackOption`.
@@ -62,27 +86,14 @@ class BlackScholes(Module):
         >>> m = BlackScholes(derivative)
         >>> m
         BSLookbackOption(strike=1.0300)
-        >>> m.inputs()
-        ['log_moneyness', 'max_log_moneyness', 'time_to_maturity', 'volatility']
-        >>> input = torch.tensor([
-        ...     [-0.01, -0.01, 0.1, 0.2],
-        ...     [ 0.00,  0.00, 0.1, 0.2],
-        ...     [ 0.01,  0.01, 0.1, 0.2]])
-        >>> m(input)
-        tensor([[...],
-                [...],
-                [...]])
     """
 
     inputs: Callable[..., List[str]]  # inputs(self) -> List[str]
     price: Callable[..., Tensor]  # price(self, ...) -> Tensor
     delta: Callable[..., Tensor]  # delta(self, ...) -> Tensor
     gamma: Callable[..., Tensor]  # gamma(self, ...) -> Tensor
+    vega: Callable[..., Tensor]  # vega(self, ...) -> Tensor
+    theta: Callable[..., Tensor]  # theta(self, ...) -> Tensor
 
     def __new__(cls, derivative):
-        return {
-            "EuropeanOption": BSEuropeanOption,
-            "LookbackOption": BSLookbackOption,
-            "AmericanBinaryOption": BSAmericanBinaryOption,
-            "EuropeanBinaryOption": BSEuropeanBinaryOption,
-        }[derivative.__class__.__name__].from_derivative(derivative)
+        return BlackScholesModuleFactory().get_class_from_derivative(derivative)
diff --git a/pfhedge/nn/modules/bs/european.py b/pfhedge/nn/modules/bs/european.py
index b1e510db..69720596 100644
--- a/pfhedge/nn/modules/bs/european.py
+++ b/pfhedge/nn/modules/bs/european.py
@@ -1,16 +1,22 @@
+from typing import Optional
+
 import torch
 from torch import Tensor
-from torch.distributions.utils import broadcast_all
 
 from pfhedge._utils.bisect import find_implied_volatility
 from pfhedge._utils.doc import _set_attr_and_docstring
 from pfhedge._utils.str import _format_float
-from pfhedge.nn.functional import d1
-from pfhedge.nn.functional import d2
-from pfhedge.nn.functional import ncdf
-from pfhedge.nn.functional import npdf
+from pfhedge.instruments import EuropeanOption
+from pfhedge.nn.functional import bs_european_delta
+from pfhedge.nn.functional import bs_european_gamma
+from pfhedge.nn.functional import bs_european_price
+from pfhedge.nn.functional import bs_european_theta
+from pfhedge.nn.functional import bs_european_vega
 
 from ._base import BSModuleMixin
+from ._base import acquire_params_from_derivative_0
+from ._base import acquire_params_from_derivative_1
+from .black_scholes import BlackScholesModuleFactory
 
 
 class BSEuropeanOption(BSModuleMixin):
@@ -55,10 +61,16 @@ class BSEuropeanOption(BSModuleMixin):
                 [0.5752]])
     """
 
-    def __init__(self, call: bool = True, strike: float = 1.0):
+    def __init__(
+        self,
+        call: bool = True,
+        strike: float = 1.0,
+        derivative: Optional[EuropeanOption] = None,
+    ) -> None:
         super().__init__()
         self.call = call
         self.strike = strike
+        self.derivative = derivative
 
     @classmethod
     def from_derivative(cls, derivative):
@@ -80,7 +92,9 @@ def from_derivative(cls, derivative):
             >>> m
             BSEuropeanOption(call=False, strike=1.)
         """
-        return cls(call=derivative.call, strike=derivative.strike)
+        return cls(
+            call=derivative.call, strike=derivative.strike, derivative=derivative
+        )
 
     def extra_repr(self) -> str:
         params = []
@@ -90,14 +104,17 @@ def extra_repr(self) -> str:
         return ", ".join(params)
 
     def delta(
-        self, log_moneyness: Tensor, time_to_maturity: Tensor, volatility: Tensor
+        self,
+        log_moneyness: Optional[Tensor] = None,
+        time_to_maturity: Optional[Tensor] = None,
+        volatility: Optional[Tensor] = None,
     ) -> Tensor:
         """Returns delta of the derivative.
 
         Args:
-            log_moneyness (torch.Tensor): Log moneyness of the underlying asset.
-            time_to_maturity (torch.Tensor): Time to expiry of the option.
-            volatility (torch.Tensor): Volatility of the underlying asset.
+            log_moneyness (torch.Tensor, optional): Log moneyness of the underlying asset.
+            time_to_maturity (torch.Tensor, optional): Time to expiry of the option.
+            volatility (torch.Tensor, optional): Volatility of the underlying asset.
 
         Shape:
             - log_moneyness: :math:`(N, *)` where
@@ -108,23 +125,36 @@ def delta(
 
         Returns:
             torch.Tensor
-        """
-        s, t, v = broadcast_all(log_moneyness, time_to_maturity, volatility)
-
-        delta = ncdf(d1(s, t, v))
-        delta = delta - 1 if not self.call else delta  # put-call parity
 
-        return delta
+        Note:
+            Parameters are not optional if the module has not accepted a derivative in its initialization.
+        """
+        (
+            log_moneyness,
+            time_to_maturity,
+            volatility,
+        ) = acquire_params_from_derivative_1(
+            self.derivative, log_moneyness, time_to_maturity, volatility
+        )
+        return bs_european_delta(
+            log_moneyness=log_moneyness,
+            time_to_maturity=time_to_maturity,
+            volatility=volatility,
+            call=self.call,
+        )
 
     def gamma(
-        self, log_moneyness: Tensor, time_to_maturity: Tensor, volatility: Tensor
+        self,
+        log_moneyness: Optional[Tensor] = None,
+        time_to_maturity: Optional[Tensor] = None,
+        volatility: Optional[Tensor] = None,
     ) -> Tensor:
         """Returns gamma of the derivative.
 
         Args:
-            log_moneyness: (torch.Tensor): Log moneyness of the underlying asset.
-            time_to_maturity (torch.Tensor): Time to expiry of the option.
-            volatility (torch.Tensor): Volatility of the underlying asset.
+            log_moneyness (torch.Tensor, optional): Log moneyness of the underlying asset.
+            time_to_maturity (torch.Tensor, optional): Time to expiry of the option.
+            volatility (torch.Tensor, optional): Volatility of the underlying asset.
 
         Shape:
             - log_moneyness: :math:`(N, *)` where
@@ -135,28 +165,36 @@ def gamma(
 
         Returns:
             torch.Tensor
-        """
 
-        s, t, v = broadcast_all(log_moneyness, time_to_maturity, volatility)
-        price = self.strike * s.exp()
-        numerator = npdf(d1(s, t, v))
-        denominator = price * v * t.sqrt()
-        output = numerator / denominator
-        return torch.where(
-            (numerator == 0).logical_and(denominator == 0),
-            torch.zeros_like(output),
-            output,
+        Note:
+            args are not optional if it doesn't accept derivative in this initialization.
+        """
+        (
+            log_moneyness,
+            time_to_maturity,
+            volatility,
+        ) = acquire_params_from_derivative_1(
+            self.derivative, log_moneyness, time_to_maturity, volatility
+        )
+        return bs_european_gamma(
+            log_moneyness=log_moneyness,
+            time_to_maturity=time_to_maturity,
+            volatility=volatility,
+            strike=self.strike,
         )
 
     def vega(
-        self, log_moneyness: Tensor, time_to_maturity: Tensor, volatility: Tensor
+        self,
+        log_moneyness: Optional[Tensor] = None,
+        time_to_maturity: Optional[Tensor] = None,
+        volatility: Optional[Tensor] = None,
     ) -> Tensor:
         """Returns vega of the derivative.
 
         Args:
-            log_moneyness: (torch.Tensor): Log moneyness of the underlying asset.
-            time_to_maturity (torch.Tensor): Time to expiry of the option.
-            volatility (torch.Tensor): Volatility of the underlying asset.
+            log_moneyness (torch.Tensor, optional): Log moneyness of the underlying asset.
+            time_to_maturity (torch.Tensor, optional): Time to expiry of the option.
+            volatility (torch.Tensor, optional): Volatility of the underlying asset.
 
         Shape:
             - log_moneyness: :math:`(N, *)` where
@@ -167,22 +205,36 @@ def vega(
 
         Returns:
             torch.Tensor
-        """
-        s, t, v = broadcast_all(log_moneyness, time_to_maturity, volatility)
-        price = self.strike * s.exp()
-        vega = npdf(d1(s, t, v)) * price * t.sqrt()
 
-        return vega
+        Note:
+            args are not optional if it doesn't accept derivative in this initialization.
+        """
+        (
+            log_moneyness,
+            time_to_maturity,
+            volatility,
+        ) = acquire_params_from_derivative_1(
+            self.derivative, log_moneyness, time_to_maturity, volatility
+        )
+        return bs_european_vega(
+            log_moneyness=log_moneyness,
+            time_to_maturity=time_to_maturity,
+            volatility=volatility,
+            strike=self.strike,
+        )
 
     def theta(
-        self, log_moneyness: Tensor, time_to_maturity: Tensor, volatility: Tensor
+        self,
+        log_moneyness: Optional[Tensor] = None,
+        time_to_maturity: Optional[Tensor] = None,
+        volatility: Optional[Tensor] = None,
     ) -> Tensor:
         """Returns theta of the derivative.
 
         Args:
-            log_moneyness: (torch.Tensor): Log moneyness of the underlying asset.
-            time_to_maturity (torch.Tensor): Time to expiry of the option.
-            volatility (torch.Tensor): Volatility of the underlying asset.
+            log_moneyness (torch.Tensor, optional): Log moneyness of the underlying asset.
+            time_to_maturity (torch.Tensor, optional): Time to expiry of the option.
+            volatility (torch.Tensor, optional): Volatility of the underlying asset.
 
         Shape:
             - log_moneyness: :math:`(N, *)` where
@@ -196,27 +248,36 @@ def theta(
 
         Returns:
             torch.Tensor
+
+        Note:
+            args are not optional if it doesn't accept derivative in this initialization.
         """
-        s, t, v = broadcast_all(log_moneyness, time_to_maturity, volatility)
-        price = self.strike * s.exp()
-        numerator = -npdf(d1(s, t, v)) * price * v
-        denominator = 2 * t.sqrt()
-        output = numerator / denominator
-        return torch.where(
-            (numerator == 0).logical_and(denominator == 0),
-            torch.zeros_like(output),
-            output,
+        (
+            log_moneyness,
+            time_to_maturity,
+            volatility,
+        ) = acquire_params_from_derivative_1(
+            self.derivative, log_moneyness, time_to_maturity, volatility
+        )
+        return bs_european_theta(
+            log_moneyness=log_moneyness,
+            time_to_maturity=time_to_maturity,
+            volatility=volatility,
+            strike=self.strike,
         )
 
     def price(
-        self, log_moneyness: Tensor, time_to_maturity: Tensor, volatility: Tensor
+        self,
+        log_moneyness: Optional[Tensor] = None,
+        time_to_maturity: Optional[Tensor] = None,
+        volatility: Optional[Tensor] = None,
     ) -> Tensor:
         """Returns price of the derivative.
 
         Args:
-            log_moneyness (torch.Tensor): Log moneyness of the underlying asset.
-            time_to_maturity (torch.Tensor): Time to expiry of the option.
-            volatility (torch.Tensor): Volatility of the underlying asset.
+            log_moneyness (torch.Tensor, optional): Log moneyness of the underlying asset.
+            time_to_maturity (torch.Tensor, optional): Time to expiry of the option.
+            volatility (torch.Tensor, optional): Volatility of the underlying asset.
 
         Shape:
             - log_moneyness: :math:`(N, *)` where
@@ -227,31 +288,37 @@ def price(
 
         Returns:
             torch.Tensor
-        """
-        s, t, v = broadcast_all(log_moneyness, time_to_maturity, volatility)
 
-        n1 = ncdf(d1(s, t, v))
-        n2 = ncdf(d2(s, t, v))
-
-        price = self.strike * (s.exp() * n1 - n2)
-
-        if not self.call:
-            price += self.strike * (1 - s.exp())  # put-call parity
-
-        return price
+        Note:
+            args are not optional if it doesn't accept derivative in this initialization.
+        """
+        (
+            log_moneyness,
+            time_to_maturity,
+            volatility,
+        ) = acquire_params_from_derivative_1(
+            self.derivative, log_moneyness, time_to_maturity, volatility
+        )
+        return bs_european_price(
+            log_moneyness=log_moneyness,
+            time_to_maturity=time_to_maturity,
+            volatility=volatility,
+            strike=self.strike,
+            call=self.call,
+        )
 
     def implied_volatility(
         self,
-        log_moneyness: Tensor,
-        time_to_maturity: Tensor,
-        price: Tensor,
+        log_moneyness: Optional[Tensor] = None,
+        time_to_maturity: Optional[Tensor] = None,
+        price: Optional[Tensor] = None,
         precision: float = 1e-6,
     ) -> Tensor:
         """Returns implied volatility of the derivative.
 
         Args:
-            log_moneyness (torch.Tensor): Log moneyness of the underlying asset.
-            time_to_maturity (torch.Tensor): Time to expiry of the option.
+            log_moneyness (torch.Tensor, optional): Log moneyness of the underlying asset.
+            time_to_maturity (torch.Tensor, optional): Time to expiry of the option.
             price (torch.Tensor): Price of the derivative.
             precision (float): Computational precision of the
                 implied volatility.
@@ -265,7 +332,19 @@ def implied_volatility(
 
         Returns
             torch.Tensor
+
+        Note:
+            args are not optional if it doesn't accept derivative in this initialization.
+            price seems optional in typing, but it isn't. It is set for the compatibility to the previous versions.
         """
+        (log_moneyness, time_to_maturity) = acquire_params_from_derivative_0(
+            self.derivative, log_moneyness, time_to_maturity
+        )
+        if price is None:
+            raise ValueError(
+                "price is required in this method. None is set only for compatibility to the previous versions."
+            )
+
         return find_implied_volatility(
             self.price,
             price=price,
@@ -275,6 +354,9 @@ def implied_volatility(
         )
 
 
+factory = BlackScholesModuleFactory()
+factory.register_module("EuropeanOption", BSEuropeanOption)
+
 # Assign docstrings so they appear in Sphinx documentation
 _set_attr_and_docstring(BSEuropeanOption, "inputs", BSModuleMixin.inputs)
 _set_attr_and_docstring(BSEuropeanOption, "forward", BSModuleMixin.forward)
diff --git a/pfhedge/nn/modules/bs/european_binary.py b/pfhedge/nn/modules/bs/european_binary.py
index 60914f25..949b8a56 100644
--- a/pfhedge/nn/modules/bs/european_binary.py
+++ b/pfhedge/nn/modules/bs/european_binary.py
@@ -1,19 +1,21 @@
 from typing import List
+from typing import Optional
 
 import torch
 from torch import Tensor
-from torch.distributions.utils import broadcast_all
 
 from pfhedge._utils.bisect import find_implied_volatility
 from pfhedge._utils.doc import _set_attr_and_docstring
 from pfhedge._utils.doc import _set_docstring
 from pfhedge._utils.str import _format_float
-from pfhedge.nn.functional import d1
-from pfhedge.nn.functional import d2
-from pfhedge.nn.functional import ncdf
-from pfhedge.nn.functional import npdf
+from pfhedge.instruments import EuropeanBinaryOption
+from pfhedge.nn.functional import bs_european_binary_delta
+from pfhedge.nn.functional import bs_european_binary_price
 
 from ._base import BSModuleMixin
+from ._base import acquire_params_from_derivative_0
+from ._base import acquire_params_from_derivative_1
+from .black_scholes import BlackScholesModuleFactory
 
 
 class BSEuropeanBinaryOption(BSModuleMixin):
@@ -58,10 +60,16 @@ class BSEuropeanBinaryOption(BSModuleMixin):
                 [6.1953]])
     """
 
-    def __init__(self, call: bool = True, strike: float = 1.0):
+    def __init__(
+        self,
+        call: bool = True,
+        strike: float = 1.0,
+        derivative: Optional[EuropeanBinaryOption] = None,
+    ):
         super().__init__()
         self.call = call
         self.strike = strike
+        self.derivative = derivative
 
     @classmethod
     def from_derivative(cls, derivative):
@@ -83,7 +91,9 @@ def from_derivative(cls, derivative):
             >>> m
             BSEuropeanBinaryOption(strike=1.1000)
         """
-        return cls(call=derivative.call, strike=derivative.strike)
+        return cls(
+            call=derivative.call, strike=derivative.strike, derivative=derivative
+        )
 
     def extra_repr(self) -> str:
         params = []
@@ -96,14 +106,17 @@ def inputs(self) -> List[str]:
         return ["log_moneyness", "time_to_maturity", "volatility"]
 
     def price(
-        self, log_moneyness: Tensor, time_to_maturity: Tensor, volatility: Tensor
+        self,
+        log_moneyness: Optional[Tensor] = None,
+        time_to_maturity: Optional[Tensor] = None,
+        volatility: Optional[Tensor] = None,
     ) -> Tensor:
         """Returns price of the derivative.
 
         Args:
-            log_moneyness (torch.Tensor): Log moneyness of the underlying asset.
-            time_to_maturity (torch.Tensor): Time to expiry of the option.
-            volatility (torch.Tensor): Volatility of the underlying asset.
+            log_moneyness (torch.Tensor, optional)): Log moneyness of the underlying asset.
+            time_to_maturity (torch.Tensor, optional)): Time to expiry of the option.
+            volatility (torch.Tensor, optional)): Volatility of the underlying asset.
 
         Shape:
             - log_moneyness: :math:`(N, *)`
@@ -113,24 +126,37 @@ def price(
 
         Returns:
             torch.Tensor
-        """
-        s, t, v = broadcast_all(log_moneyness, time_to_maturity, volatility)
-
-        price = ncdf(d2(s, t, v))
-        price = 1.0 - price if not self.call else price  # put-call parity
 
-        return price
+        Note:
+            args are not optional if it doesn't accept derivative in this initialization.
+        """
+        (
+            log_moneyness,
+            time_to_maturity,
+            volatility,
+        ) = acquire_params_from_derivative_1(
+            self.derivative, log_moneyness, time_to_maturity, volatility
+        )
+        return bs_european_binary_price(
+            log_moneyness=log_moneyness,
+            time_to_maturity=time_to_maturity,
+            volatility=volatility,
+            call=self.call,
+        )
 
     @torch.enable_grad()
     def delta(
-        self, log_moneyness: Tensor, time_to_maturity: Tensor, volatility: Tensor
+        self,
+        log_moneyness: Optional[Tensor] = None,
+        time_to_maturity: Optional[Tensor] = None,
+        volatility: Optional[Tensor] = None,
     ) -> Tensor:
         """Returns delta of the derivative.
 
         Args:
-            log_moneyness: (torch.Tensor): Log moneyness of the underlying asset.
-            time_to_maturity (torch.Tensor): Time to expiry of the option.
-            volatility (torch.Tensor): Volatility of the underlying asset.
+            log_moneyness: (torch.Tensor, optional): Log moneyness of the underlying asset.
+            time_to_maturity (torch.Tensor, optional): Time to expiry of the option.
+            volatility (torch.Tensor, optional): Volatility of the underlying asset.
 
         Shape:
             - log_moneyness: :math:`(N, *)`
@@ -140,32 +166,37 @@ def delta(
 
         Returns:
             torch.Tensor
-        """
-        s, t, v = broadcast_all(log_moneyness, time_to_maturity, volatility)
 
-        spot = s.exp() * self.strike
-
-        numerator = npdf(d2(s, t, v))
-        denominator = spot * v * t.sqrt()
-        delta = numerator / denominator
-        delta = torch.where(
-            (numerator == 0).logical_and(denominator == 0),
-            torch.zeros_like(delta),
-            delta,
+        Note:
+            Parameters are not optional if the module has not accepted a derivative in its initialization.
+        """
+        (
+            log_moneyness,
+            time_to_maturity,
+            volatility,
+        ) = acquire_params_from_derivative_1(
+            self.derivative, log_moneyness, time_to_maturity, volatility
+        )
+        return bs_european_binary_delta(
+            log_moneyness=log_moneyness,
+            time_to_maturity=time_to_maturity,
+            volatility=volatility,
+            call=self.call,
+            strike=self.strike,
         )
-        delta = -delta if not self.call else delta  # put-call parity
-
-        return delta
 
     def gamma(
-        self, log_moneyness: Tensor, time_to_maturity: Tensor, volatility: Tensor
+        self,
+        log_moneyness: Optional[Tensor] = None,
+        time_to_maturity: Optional[Tensor] = None,
+        volatility: Optional[Tensor] = None,
     ) -> Tensor:
         """Returns gamma of the derivative.
 
         Args:
-            log_moneyness (torch.Tensor): Log moneyness of the underlying asset.
-            time_to_maturity (torch.Tensor): Time to expiry of the option.
-            volatility (torch.Tensor): Volatility of the underlying asset.
+            log_moneyness (torch.Tensor, optional): Log moneyness of the underlying asset.
+            time_to_maturity (torch.Tensor, optional): Time to expiry of the option.
+            volatility (torch.Tensor, optional): Volatility of the underlying asset.
 
         Shape:
             - log_moneyness: :math:`(N, *)`
@@ -175,24 +206,37 @@ def gamma(
 
         Returns:
             torch.Tensor
+
+        Note:
+            args are not optional if it doesn't accept derivative in this initialization.
         """
+        (
+            log_moneyness,
+            time_to_maturity,
+            volatility,
+        ) = acquire_params_from_derivative_1(
+            self.derivative, log_moneyness, time_to_maturity, volatility
+        )
         # TODO(simaki): Directly compute gamma.
         return super().gamma(
-            strike=self.strike,
             log_moneyness=log_moneyness,
             time_to_maturity=time_to_maturity,
             volatility=volatility,
+            strike=self.strike,
         )
 
     def vega(
-        self, log_moneyness: Tensor, time_to_maturity: Tensor, volatility: Tensor
+        self,
+        log_moneyness: Optional[Tensor] = None,
+        time_to_maturity: Optional[Tensor] = None,
+        volatility: Optional[Tensor] = None,
     ) -> Tensor:
         """Returns vega of the derivative.
 
         Args:
-            log_moneyness (torch.Tensor): Log moneyness of the underlying asset.
-            time_to_maturity (torch.Tensor): Time to expiry of the option.
-            volatility (torch.Tensor): Volatility of the underlying asset.
+            log_moneyness (torch.Tensor, optional): Log moneyness of the underlying asset.
+            time_to_maturity (torch.Tensor, optional): Time to expiry of the option.
+            volatility (torch.Tensor, optional): Volatility of the underlying asset.
 
         Shape:
             - log_moneyness: :math:`(N, *)`
@@ -202,24 +246,37 @@ def vega(
 
         Returns:
             torch.Tensor
+
+        Note:
+            args are not optional if it doesn't accept derivative in this initialization.
         """
+        (
+            log_moneyness,
+            time_to_maturity,
+            volatility,
+        ) = acquire_params_from_derivative_1(
+            self.derivative, log_moneyness, time_to_maturity, volatility
+        )
         # TODO: Directly compute theta.
         return super().vega(
-            strike=self.strike,
             log_moneyness=log_moneyness,
             time_to_maturity=time_to_maturity,
             volatility=volatility,
+            strike=self.strike,
         )
 
     def theta(
-        self, log_moneyness: Tensor, time_to_maturity: Tensor, volatility: Tensor
+        self,
+        log_moneyness: Optional[Tensor] = None,
+        time_to_maturity: Optional[Tensor] = None,
+        volatility: Optional[Tensor] = None,
     ) -> Tensor:
         """Returns theta of the derivative.
 
         Args:
-            log_moneyness (torch.Tensor): Log moneyness of the underlying asset.
-            time_to_maturity (torch.Tensor): Time to expiry of the option.
-            volatility (torch.Tensor): Volatility of the underlying asset.
+            log_moneyness (torch.Tensor, optional): Log moneyness of the underlying asset.
+            time_to_maturity (torch.Tensor, optional): Time to expiry of the option.
+            volatility (torch.Tensor, optional): Volatility of the underlying asset.
 
         Shape:
             - log_moneyness: :math:`(N, *)`
@@ -232,37 +289,60 @@ def theta(
 
         Returns:
             torch.Tensor
+
+        Note:
+            args are not optional if it doesn't accept derivative in this initialization.
         """
+        (
+            log_moneyness,
+            time_to_maturity,
+            volatility,
+        ) = acquire_params_from_derivative_1(
+            self.derivative, log_moneyness, time_to_maturity, volatility
+        )
         # TODO: Directly compute theta.
         return super().theta(
-            strike=self.strike,
             log_moneyness=log_moneyness,
             time_to_maturity=time_to_maturity,
             volatility=volatility,
+            strike=self.strike,
         )
 
     def implied_volatility(
         self,
-        log_moneyness: Tensor,
-        time_to_maturity: Tensor,
-        price: Tensor,
+        log_moneyness: Optional[Tensor] = None,
+        time_to_maturity: Optional[Tensor] = None,
+        price: Optional[Tensor] = None,
         precision: float = 1e-6,
     ) -> Tensor:
         """Returns implied volatility of the derivative.
 
         Args:
-            log_moneyness (torch.Tensor): Log moneyness of the underlying asset.
-            time_to_maturity (torch.Tensor): Time to expiry of the option.
+            log_moneyness (torch.Tensor, optional): Log moneyness of the underlying asset.
+            time_to_maturity (torch.Tensor, optional): Time to expiry of the option.
             price (torch.Tensor): Price of the derivative.
+            precision (float): Computational precision of the
+                implied volatility.
 
         Shape:
             - log_moneyness: :math:`(N, *)`
             - time_to_maturity: :math:`(N, *)`
+            - price: :math:`(N, *)`
             - output: :math:`(N, *)`
 
         Returns:
             torch.Tensor
+        Note:
+            args are not optional if it doesn't accept derivative in this initialization.
+            price seems optional in typing, but it isn't. It is set for the compatibility to the previous versions.
         """
+        (log_moneyness, time_to_maturity) = acquire_params_from_derivative_0(
+            self.derivative, log_moneyness, time_to_maturity
+        )
+        if price is None:
+            raise ValueError(
+                "price is required in this method. None is set only for compatibility to the previous versions."
+            )
         return find_implied_volatility(
             self.price,
             price=price,
@@ -272,6 +352,9 @@ def implied_volatility(
         )
 
 
+factory = BlackScholesModuleFactory()
+factory.register_module("EuropeanBinaryOption", BSEuropeanBinaryOption)
+
 # Assign docstrings so they appear in Sphinx documentation
 _set_docstring(BSEuropeanBinaryOption, "inputs", BSModuleMixin.inputs)
 _set_attr_and_docstring(BSEuropeanBinaryOption, "forward", BSModuleMixin.forward)
diff --git a/pfhedge/nn/modules/bs/lookback.py b/pfhedge/nn/modules/bs/lookback.py
index 72d83551..3c868a11 100644
--- a/pfhedge/nn/modules/bs/lookback.py
+++ b/pfhedge/nn/modules/bs/lookback.py
@@ -1,4 +1,5 @@
 from typing import List
+from typing import Optional
 
 import torch
 from torch import Tensor
@@ -7,19 +8,27 @@
 from pfhedge._utils.doc import _set_attr_and_docstring
 from pfhedge._utils.doc import _set_docstring
 from pfhedge._utils.str import _format_float
-from pfhedge.nn.functional import d1 as compute_d1
-from pfhedge.nn.functional import d2 as compute_d2
-from pfhedge.nn.functional import ncdf
-from pfhedge.nn.functional import npdf
+from pfhedge.instruments import LookbackOption
+from pfhedge.nn.functional import bs_lookback_delta
+from pfhedge.nn.functional import bs_lookback_gamma
+from pfhedge.nn.functional import bs_lookback_price
+from pfhedge.nn.functional import bs_lookback_theta
+from pfhedge.nn.functional import bs_lookback_vega
 
 from ._base import BSModuleMixin
+from ._base import acquire_params_from_derivative_0
+from ._base import acquire_params_from_derivative_2
+from .black_scholes import BlackScholesModuleFactory
 
 
 class BSLookbackOption(BSModuleMixin):
     """Black-Scholes formula for a lookback option with a fixed strike.
 
     Note:
-        Risk-free rate is set to zero.
+        - The formulas are for continuous monitoring while
+          :class:`pfhedge.instruments.LookbackOption` monitors spot prices discretely.
+          To get adjustment for discrete monitoring, see, for instance,
+          Broadie, Glasserman, and Kou (1999).
 
     .. seealso::
         - :class:`pfhedge.nn.BlackScholes`:
@@ -30,6 +39,9 @@ class BSLookbackOption(BSModuleMixin):
     References:
         - Conze, A., 1991. Path dependent options: The case of lookback options.
           The Journal of Finance, 46(5), pp.1893-1907.
+        - Broadie, M., Glasserman, P. and Kou, S.G., 1999.
+          Connecting discrete and continuous path-dependent options.
+          Finance and Stochastics, 3(1), pp.55-82.
 
     Args:
         call (bool, default=True): Specifies whether the option is call or put.
@@ -58,7 +70,12 @@ class BSLookbackOption(BSModuleMixin):
                 [1.0515]])
     """
 
-    def __init__(self, call: bool = True, strike: float = 1.0) -> None:
+    def __init__(
+        self,
+        call: bool = True,
+        strike: float = 1.0,
+        derivative: Optional[LookbackOption] = None,
+    ) -> None:
         if not call:
             raise ValueError(
                 f"{self.__class__.__name__} for a put option is not yet supported."
@@ -67,6 +84,7 @@ def __init__(self, call: bool = True, strike: float = 1.0) -> None:
         super().__init__()
         self.call = call
         self.strike = strike
+        self.derivative = derivative
 
     @classmethod
     def from_derivative(cls, derivative):
@@ -88,7 +106,9 @@ def from_derivative(cls, derivative):
             >>> m
             BSLookbackOption(strike=1.1000)
         """
-        return cls(call=derivative.call, strike=derivative.strike)
+        return cls(
+            call=derivative.call, strike=derivative.strike, derivative=derivative
+        )
 
     def extra_repr(self) -> str:
         params = []
@@ -100,10 +120,10 @@ def inputs(self) -> List[str]:
 
     def price(
         self,
-        log_moneyness: Tensor,
-        max_log_moneyness: Tensor,
-        time_to_maturity: Tensor,
-        volatility: Tensor,
+        log_moneyness: Optional[Tensor] = None,
+        max_log_moneyness: Optional[Tensor] = None,
+        time_to_maturity: Optional[Tensor] = None,
+        volatility: Optional[Tensor] = None,
     ) -> Tensor:
         r"""Returns price of the derivative.
 
@@ -126,10 +146,10 @@ def price(
         :math:`d_2' = [\log(S(0) / M) - \frac12 \sigma^2 T] / \sigma \sqrt{T}`.
 
         Args:
-            log_moneyness (torch.Tensor): Log moneyness of the underlying asset.
-            max_log_moneyness (torch.Tensor): Cumulative maximum of the log moneyness.
-            time_to_maturity (torch.Tensor): Time to expiry of the option.
-            volatility (torch.Tensor): Volatility of the underlying asset.
+            log_moneyness (torch.Tensor, optional): Log moneyness of the underlying asset.
+            max_log_moneyness (torch.Tensor, optional): Cumulative maximum of the log moneyness.
+            time_to_maturity (torch.Tensor, optional): Time to expiry of the option.
+            volatility (torch.Tensor, optional): Volatility of the underlying asset.
 
         Shape:
             - log_moneyness: :math:`(N, *)` where
@@ -141,48 +161,46 @@ def price(
 
         Returns:
             torch.Tensor
+
+        Note:
+            Parameters are not optional if the module has not accepted a derivative in its initialization.
         """
-        s, m, t, v = map(
-            torch.as_tensor,
-            (log_moneyness, max_log_moneyness, time_to_maturity, volatility),
+        (
+            log_moneyness,
+            max_log_moneyness,
+            time_to_maturity,
+            volatility,
+        ) = acquire_params_from_derivative_2(
+            self.derivative,
+            log_moneyness,
+            max_log_moneyness,
+            time_to_maturity,
+            volatility,
         )
-
-        spot = s.exp() * self.strike
-        max = m.exp() * self.strike
-        d1 = compute_d1(s, t, v)
-        d2 = compute_d2(s, t, v)
-        m1 = compute_d1(s - m, t, v)  # d' in the paper
-        m2 = compute_d2(s - m, t, v)
-
-        # when max < strike
-        price_0 = spot * (
-            ncdf(d1) + v * t.sqrt() * (d1 * ncdf(d1) + npdf(d1))
-        ) - self.strike * ncdf(d2)
-        # when max >= strike
-        price_1 = (
-            spot * (ncdf(m1) + v * t.sqrt() * (m1 * ncdf(m1) + npdf(m1)))
-            - self.strike
-            + max * (1 - ncdf(m2))
+        return bs_lookback_price(
+            log_moneyness=log_moneyness,
+            max_log_moneyness=max_log_moneyness,
+            time_to_maturity=time_to_maturity,
+            volatility=volatility,
+            strike=self.strike,
         )
 
-        return torch.where(max < self.strike, price_0, price_1)
-
     @torch.enable_grad()
     def delta(
         self,
-        log_moneyness: Tensor,
-        max_log_moneyness: Tensor,
-        time_to_maturity: Tensor,
-        volatility: Tensor,
+        log_moneyness: Optional[Tensor] = None,
+        max_log_moneyness: Optional[Tensor] = None,
+        time_to_maturity: Optional[Tensor] = None,
+        volatility: Optional[Tensor] = None,
         create_graph: bool = False,
     ) -> Tensor:
         """Returns delta of the derivative.
 
         Args:
-            log_moneyness (torch.Tensor): Log moneyness of the underlying asset.
-            max_log_moneyness (torch.Tensor): Cumulative maximum of the log moneyness.
-            time_to_maturity (torch.Tensor): Time to expiry of the option.
-            volatility (torch.Tensor): Volatility of the underlying asset.
+            log_moneyness (torch.Tensor, optional): Log moneyness of the underlying asset.
+            max_log_moneyness (torch.Tensor, optional): Cumulative maximum of the log moneyness.
+            time_to_maturity (torch.Tensor, optional): Time to expiry of the option.
+            volatility (torch.Tensor, optional): Volatility of the underlying asset.
             create_graph (bool, default=False): If True, graph of the derivative
                 will be constructed. This option is used to compute gamma.
 
@@ -196,30 +214,46 @@ def delta(
 
         Returns:
             torch.Tensor
+
+        Note:
+            Parameters are not optional if the module has not accepted a derivative in its initialization.
         """
+        (
+            log_moneyness,
+            max_log_moneyness,
+            time_to_maturity,
+            volatility,
+        ) = acquire_params_from_derivative_2(
+            self.derivative,
+            log_moneyness,
+            max_log_moneyness,
+            time_to_maturity,
+            volatility,
+        )
         return super().delta(
             log_moneyness=log_moneyness,
             max_log_moneyness=max_log_moneyness,
             time_to_maturity=time_to_maturity,
             volatility=volatility,
             create_graph=create_graph,
+            strike=self.strike,
         )
 
     @torch.enable_grad()
     def gamma(
         self,
-        log_moneyness: Tensor,
-        max_log_moneyness: Tensor,
-        time_to_maturity: Tensor,
-        volatility: Tensor,
+        log_moneyness: Optional[Tensor] = None,
+        max_log_moneyness: Optional[Tensor] = None,
+        time_to_maturity: Optional[Tensor] = None,
+        volatility: Optional[Tensor] = None,
     ) -> Tensor:
         """Returns gamma of the derivative.
 
         Args:
-            log_moneyness (torch.Tensor): Log moneyness of the underlying asset.
-            max_log_moneyness (torch.Tensor): Cumulative maximum of the log moneyness.
-            time_to_maturity (torch.Tensor): Time to expiry of the option.
-            volatility (torch.Tensor):
+            log_moneyness (torch.Tensor, optional): Log moneyness of the underlying asset.
+            max_log_moneyness (torch.Tensor, optional): Cumulative maximum of the log moneyness.
+            time_to_maturity (torch.Tensor, optional): Time to expiry of the option.
+            volatility (torch.Tensor, optional):
                 Volatility of the underlying asset.
 
         Shape:
@@ -232,7 +266,22 @@ def gamma(
 
         Returns:
             torch.Tensor
+
+        Note:
+            args are not optional if it doesn't accept derivative in this initialization.
         """
+        (
+            log_moneyness,
+            max_log_moneyness,
+            time_to_maturity,
+            volatility,
+        ) = acquire_params_from_derivative_2(
+            self.derivative,
+            log_moneyness,
+            max_log_moneyness,
+            time_to_maturity,
+            volatility,
+        )
         return super().gamma(
             strike=self.strike,
             log_moneyness=log_moneyness,
@@ -244,18 +293,18 @@ def gamma(
     @torch.enable_grad()
     def vega(
         self,
-        log_moneyness: Tensor,
-        max_log_moneyness: Tensor,
-        time_to_maturity: Tensor,
-        volatility: Tensor,
+        log_moneyness: Optional[Tensor] = None,
+        max_log_moneyness: Optional[Tensor] = None,
+        time_to_maturity: Optional[Tensor] = None,
+        volatility: Optional[Tensor] = None,
     ) -> Tensor:
         """Returns vega of the derivative.
 
         Args:
-            log_moneyness (torch.Tensor): Log moneyness of the underlying asset.
-            max_log_moneyness (torch.Tensor): Cumulative maximum of the log moneyness.
-            time_to_maturity (torch.Tensor): Time to expiry of the option.
-            volatility (torch.Tensor):
+            log_moneyness (torch.Tensor, optional): Log moneyness of the underlying asset.
+            max_log_moneyness (torch.Tensor, optional): Cumulative maximum of the log moneyness.
+            time_to_maturity (torch.Tensor, optional): Time to expiry of the option.
+            volatility (torch.Tensor, optional):
                 Volatility of the underlying asset.
 
         Shape:
@@ -268,7 +317,22 @@ def vega(
 
         Returns:
             torch.Tensor
+
+        Note:
+            args are not optional if it doesn't accept derivative in this initialization.
         """
+        (
+            log_moneyness,
+            max_log_moneyness,
+            time_to_maturity,
+            volatility,
+        ) = acquire_params_from_derivative_2(
+            self.derivative,
+            log_moneyness,
+            max_log_moneyness,
+            time_to_maturity,
+            volatility,
+        )
         return super().vega(
             strike=self.strike,
             log_moneyness=log_moneyness,
@@ -280,18 +344,18 @@ def vega(
     @torch.enable_grad()
     def theta(
         self,
-        log_moneyness: Tensor,
-        max_log_moneyness: Tensor,
-        time_to_maturity: Tensor,
-        volatility: Tensor,
+        log_moneyness: Optional[Tensor] = None,
+        max_log_moneyness: Optional[Tensor] = None,
+        time_to_maturity: Optional[Tensor] = None,
+        volatility: Optional[Tensor] = None,
     ) -> Tensor:
         """Returns theta of the derivative.
 
         Args:
-            log_moneyness (torch.Tensor): Log moneyness of the underlying asset.
-            max_log_moneyness (torch.Tensor): Cumulative maximum of the log moneyness.
-            time_to_maturity (torch.Tensor): Time to expiry of the option.
-            volatility (torch.Tensor):
+            log_moneyness (torch.Tensor, optional): Log moneyness of the underlying asset.
+            max_log_moneyness (torch.Tensor, optional): Cumulative maximum of the log moneyness.
+            time_to_maturity (torch.Tensor, optional): Time to expiry of the option.
+            volatility (torch.Tensor, optional):
                 Volatility of the underlying asset.
 
         Shape:
@@ -307,7 +371,22 @@ def theta(
 
         Returns:
             torch.Tensor
+
+        Note:
+            args are not optional if it doesn't accept derivative in this initialization.
         """
+        (
+            log_moneyness,
+            max_log_moneyness,
+            time_to_maturity,
+            volatility,
+        ) = acquire_params_from_derivative_2(
+            self.derivative,
+            log_moneyness,
+            max_log_moneyness,
+            time_to_maturity,
+            volatility,
+        )
         return super().theta(
             strike=self.strike,
             log_moneyness=log_moneyness,
@@ -318,18 +397,18 @@ def theta(
 
     def implied_volatility(
         self,
-        log_moneyness: Tensor,
-        max_log_moneyness: Tensor,
-        time_to_maturity: Tensor,
-        price: Tensor,
+        log_moneyness: Optional[Tensor] = None,
+        max_log_moneyness: Optional[Tensor] = None,
+        time_to_maturity: Optional[Tensor] = None,
+        price: Optional[Tensor] = None,
         precision: float = 1e-6,
     ) -> Tensor:
         """Returns implied volatility of the derivative.
 
         Args:
-            log_moneyness (torch.Tensor): Log moneyness of the underlying asset.
-            max_log_moneyness (torch.Tensor): Cumulative maximum of the log moneyness.
-            time_to_maturity (torch.Tensor): Time to expiry of the option.
+            log_moneyness (torch.Tensor, optional): Log moneyness of the underlying asset.
+            max_log_moneyness (torch.Tensor, optional): Cumulative maximum of the log moneyness.
+            time_to_maturity (torch.Tensor, optional): Time to expiry of the option.
             price (torch.Tensor): Price of the derivative.
             precision (float, default=1e-6): Precision of the implied volatility.
 
@@ -343,7 +422,18 @@ def implied_volatility(
 
         Returns:
             torch.Tensor
+
+        Note:
+            args are not optional if it doesn't accept derivative in this initialization.
+            price seems optional in typing, but it isn't. It is set for the compatibility to the previous versions.
         """
+        (log_moneyness, time_to_maturity) = acquire_params_from_derivative_0(
+            self.derivative, log_moneyness, time_to_maturity
+        )
+        if price is None:
+            raise ValueError(
+                "price is required in this method. None is set only for compatibility to the previous versions."
+            )
         return find_implied_volatility(
             self.price,
             price=price,
@@ -354,6 +444,9 @@ def implied_volatility(
         )
 
 
+factory = BlackScholesModuleFactory()
+factory.register_module("LookbackOption", BSLookbackOption)
+
 # Assign docstrings so they appear in Sphinx documentation
 _set_docstring(BSLookbackOption, "inputs", BSModuleMixin.inputs)
 _set_attr_and_docstring(BSLookbackOption, "forward", BSModuleMixin.forward)
diff --git a/pfhedge/nn/modules/hedger.py b/pfhedge/nn/modules/hedger.py
index c11c1884..3927cc2c 100644
--- a/pfhedge/nn/modules/hedger.py
+++ b/pfhedge/nn/modules/hedger.py
@@ -1,3 +1,4 @@
+from typing import Any
 from typing import Callable
 from typing import List
 from typing import Optional
@@ -539,7 +540,7 @@ def fit(
         """
         optimizer = self._configure_optimizer(derivative, optimizer)
 
-        def compute_loss(**kwargs) -> Tensor:
+        def compute_loss(**kwargs: Any) -> Tensor:
             return self.compute_loss(
                 derivative,
                 hedge=hedge,
diff --git a/pfhedge/stochastic/brownian.py b/pfhedge/stochastic/brownian.py
index 21715897..34a88650 100644
--- a/pfhedge/stochastic/brownian.py
+++ b/pfhedge/stochastic/brownian.py
@@ -82,7 +82,7 @@ def generate_brownian(
 
     init_value = init_state[0]
     # randn = torch.randn((n_paths, n_steps), dtype=dtype, device=device)
-    randn = engine((n_paths, n_steps), dtype=dtype, device=device)
+    randn = engine(*(n_paths, n_steps), dtype=dtype, device=device)
     randn[:, 0] = 0.0
     drift = mu * dt * torch.arange(n_steps).to(randn)
     brown = randn.new_tensor(dt).sqrt() * randn.cumsum(1)
diff --git a/pfhedge/stochastic/cir.py b/pfhedge/stochastic/cir.py
index d835f816..7340fbe2 100644
--- a/pfhedge/stochastic/cir.py
+++ b/pfhedge/stochastic/cir.py
@@ -8,7 +8,7 @@
 from pfhedge._utils.typing import TensorOrScalar
 
 
-def _get_epsilon(dtype: Optional[torch.dtype]):
+def _get_epsilon(dtype: Optional[torch.dtype]) -> float:
     return torch.finfo(dtype).tiny if dtype else torch.finfo().tiny
 
 
diff --git a/pfhedge/stochastic/engine.py b/pfhedge/stochastic/engine.py
index 112f2a9f..68db2a72 100644
--- a/pfhedge/stochastic/engine.py
+++ b/pfhedge/stochastic/engine.py
@@ -15,7 +15,7 @@ def _box_muller(input0: Tensor, input1: Tensor) -> Tuple[Tensor, Tensor]:
     return z0, z1
 
 
-def _get_numel(size: Tuple[int, ...]):
+def _get_numel(size: Tuple[int, ...]) -> int:
     out = 1
     for dim in size:
         out *= dim
@@ -42,11 +42,11 @@ def __init__(self, scramble: bool = False, seed: Optional[int] = None):
 
     def __call__(
         self,
-        *size: Tuple[int, ...],
+        *size: int,
         dtype: Optional[torch.dtype] = None,
         device: Optional[torch.device] = None,
     ) -> Tensor:
-        numel = _get_numel(*size)
+        numel = _get_numel(tuple(size))
         output = self._generate_1d(numel, dtype=dtype, device=device)
         output.resize_(*size)
         return output
diff --git a/pfhedge/stochastic/random.py b/pfhedge/stochastic/random.py
index 3bee1fbc..711e7388 100644
--- a/pfhedge/stochastic/random.py
+++ b/pfhedge/stochastic/random.py
@@ -1,9 +1,18 @@
+from typing import Optional
+
 import torch
+from torch import Tensor
 
 from .engine import RandnSobolBoxMuller
 
 
-def randn_antithetic(*size, dtype=None, device=None, dim=0, shuffle=True):
+def randn_antithetic(
+    *size: int,
+    dtype: Optional[torch.dtype] = None,
+    device: Optional[torch.device] = None,
+    dim: int = 0,
+    shuffle: bool = True
+) -> Tensor:
     """Returns a tensor filled with random numbers obtained by an antithetic sampling.
 
     The output should be a normal distribution with mean 0 and variance 1
@@ -40,8 +49,8 @@ def randn_antithetic(*size, dtype=None, device=None, dim=0, shuffle=True):
     if dim != 0:
         raise ValueError("dim != 0 is not supported.")
 
-    size = list(size)
-    size_half = [-(-size[0] // 2)] + size[1:]
+    size_list = list(size)
+    size_half = [-(-size_list[0] // 2)] + size_list[1:]
     randn = torch.randn(*size_half, dtype=dtype, device=device)
 
     output = torch.cat((randn, -randn), dim=0)
@@ -49,12 +58,18 @@ def randn_antithetic(*size, dtype=None, device=None, dim=0, shuffle=True):
     if shuffle:
         output = output[torch.randperm(output.size(dim))]
 
-    output = output[: size[0]]
+    output = output[: size_list[0]]
 
     return output
 
 
-def randn_sobol_boxmuller(*size, dtype=None, device=None, scramble=True, seed=None):
+def randn_sobol_boxmuller(
+    *size: int,
+    dtype: Optional[torch.dtype] = None,
+    device: Optional[torch.device] = None,
+    scramble: bool = True,
+    seed: Optional[int] = None
+) -> Tensor:
     """Returns a tensor filled with random numbers obtained by a Sobol sequence
     applied with the Box-Muller transformation.
 
@@ -80,7 +95,7 @@ def randn_sobol_boxmuller(*size, dtype=None, device=None, scramble=True, seed=No
         >>> from pfhedge.stochastic import randn_sobol_boxmuller
         >>>
         >>> _ = torch.manual_seed(42)
-        >>> output = randn_sobol_boxmuller((4, 3))
+        >>> output = randn_sobol_boxmuller(4, 3)
         >>> output
         tensor([[ 0.0559,  0.4954, -0.8578],
                 [-0.7492, -1.0370, -0.4778],
diff --git a/pfhedge/version.py b/pfhedge/version.py
new file mode 100644
index 00000000..1317d755
--- /dev/null
+++ b/pfhedge/version.py
@@ -0,0 +1 @@
+__version__ = "0.18.0"
diff --git a/pyproject.toml b/pyproject.toml
index d0951426..5a7055a9 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -1,6 +1,6 @@
 [tool.poetry]
 name = "pfhedge"
-version = "0.17.0"
+version = "0.18.0"
 description = "Deep Hedging in PyTorch"
 authors = ["Shota Imaki <shota.imaki.0801@gmail.com>"]
 license = "MIT"
diff --git a/tests/instruments/primary/test_brownian.py b/tests/instruments/primary/test_brownian.py
index b536f4aa..ad0f829c 100644
--- a/tests/instruments/primary/test_brownian.py
+++ b/tests/instruments/primary/test_brownian.py
@@ -146,18 +146,36 @@ def test_to_dtype(self, dtype):
         assert s.dtype == torch.float64
         assert s.spot.dtype == torch.float64
 
+        s = BrownianStock()
+        s.simulate()
+        s.float64()
+        assert s.dtype == torch.float64
+        assert s.spot.dtype == torch.float64
+
         s = BrownianStock()
         s.simulate()
         s.float()
         assert s.dtype == torch.float32
         assert s.spot.dtype == torch.float32
 
+        s = BrownianStock()
+        s.simulate()
+        s.float32()
+        assert s.dtype == torch.float32
+        assert s.spot.dtype == torch.float32
+
         s = BrownianStock()
         s.simulate()
         s.half()
         assert s.dtype == torch.float16
         assert s.spot.dtype == torch.float16
 
+        s = BrownianStock()
+        s.simulate()
+        s.float16()
+        assert s.dtype == torch.float16
+        assert s.spot.dtype == torch.float16
+
         s = BrownianStock()
         s.simulate()
         s.bfloat16()
diff --git a/tests/nn/modules/bs/test_american_binary.py b/tests/nn/modules/bs/test_american_binary.py
index ad25e78c..7d275d27 100644
--- a/tests/nn/modules/bs/test_american_binary.py
+++ b/tests/nn/modules/bs/test_american_binary.py
@@ -2,7 +2,6 @@
 
 import pytest
 import torch
-from torch import Tensor
 from torch.testing import assert_close
 
 from pfhedge.features._getter import get_feature
@@ -62,6 +61,107 @@ def test_check_delta(self):
         expect = torch.tensor([0.0])
         assert_close(result, expect)
 
+    @pytest.mark.parametrize("call", [True])
+    def test_delta_2(self, call: bool):
+        m = BSAmericanBinaryOption(call=call)
+        with pytest.raises(ValueError):
+            m.delta(
+                torch.tensor(0.0),
+                torch.tensor(0.0),
+                torch.tensor(-1.0),
+                torch.tensor(0.2),
+            )
+        with pytest.raises(ValueError):
+            m.delta(
+                torch.tensor(0.0),
+                torch.tensor(0.0),
+                torch.tensor(1.0),
+                torch.tensor(-0.2),
+            )
+
+    @pytest.mark.parametrize("call", [True])
+    def test_delta_3(self, call: bool):
+        derivative = AmericanBinaryOption(BrownianStock(), call=call)
+        m = BSAmericanBinaryOption.from_derivative(derivative)
+        m2 = BSAmericanBinaryOption(call=call)
+        with pytest.raises(AttributeError):
+            m.delta()
+        with pytest.raises(AttributeError):
+            m.delta(None, torch.tensor(1), torch.tensor(2), torch.tensor(3))
+        with pytest.raises(AttributeError):
+            m.delta(torch.tensor(1), None, torch.tensor(2), torch.tensor(3))
+        with pytest.raises(AttributeError):
+            m.delta(torch.tensor(1), torch.tensor(2), None, torch.tensor(3))
+        # ToDo: #530
+        # with pytest.raises(AttributeError):
+        #     m.delta(torch.tensor(1), torch.tensor(2), torch.tensor(3), None)
+        torch.manual_seed(42)
+        derivative.simulate(n_paths=1)
+        result = m.delta()
+        expect = m2.delta(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result, expect)
+        result = m.delta(
+            None,
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result, expect)
+        result = m.delta(
+            derivative.log_moneyness(),
+            None,
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result, expect)
+        result = m.delta(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            None,
+            derivative.underlier.volatility,
+        )
+        assert_close(result, expect)
+        result = m.delta(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            None,
+        )
+        assert_close(result, expect)
+        with pytest.raises(ValueError):
+            m2.delta(
+                None,
+                derivative.max_log_moneyness(),
+                derivative.time_to_maturity(),
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.delta(
+                derivative.log_moneyness(),
+                None,
+                derivative.time_to_maturity(),
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.delta(
+                derivative.log_moneyness(),
+                derivative.max_log_moneyness(),
+                None,
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.delta(
+                derivative.log_moneyness(),
+                derivative.max_log_moneyness(),
+                derivative.time_to_maturity(),
+                None,
+            )
+
     def test_check_gamma(self):
         m = BSAmericanBinaryOption()
 
@@ -85,6 +185,106 @@ def test_check_gamma(self):
         expect = torch.tensor([0.0])
         assert_close(result, expect)
 
+    @pytest.mark.parametrize("call", [True])
+    def test_gamma_2(self, call: bool):
+        m = BSAmericanBinaryOption(call=call)
+        with pytest.raises(ValueError):
+            m.gamma(
+                torch.tensor(0.0),
+                torch.tensor(0.0),
+                torch.tensor(-1.0),
+                torch.tensor(0.2),
+            )
+        with pytest.raises(ValueError):
+            m.gamma(
+                torch.tensor(0.0),
+                torch.tensor(0.0),
+                torch.tensor(1.0),
+                torch.tensor(-0.2),
+            )
+
+    @pytest.mark.parametrize("call", [True])
+    def test_gamma_3(self, call: bool):
+        derivative = AmericanBinaryOption(BrownianStock(), call=call)
+        m = BSAmericanBinaryOption.from_derivative(derivative)
+        m2 = BSAmericanBinaryOption(call=call)
+        with pytest.raises(AttributeError):
+            m.gamma(None, torch.tensor(1), torch.tensor(2), torch.tensor(3))
+        with pytest.raises(AttributeError):
+            m.gamma(torch.tensor(1), None, torch.tensor(2), torch.tensor(3))
+        with pytest.raises(AttributeError):
+            m.gamma(torch.tensor(1), torch.tensor(2), None, torch.tensor(3))
+        # ToDo: #530
+        # with pytest.raises(AttributeError):
+        #     m.gamma(torch.tensor(1), torch.tensor(2), torch.tensor(3), None)
+        torch.manual_seed(42)
+        derivative.simulate(n_paths=1)
+        result = m.gamma()
+        expect = m2.gamma(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        # ToDo: [..., :-1] should be removed
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.gamma(
+            None,
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.gamma(
+            derivative.log_moneyness(),
+            None,
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.gamma(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            None,
+            derivative.underlier.volatility,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.gamma(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            None,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        with pytest.raises(ValueError):
+            m2.gamma(
+                None,
+                derivative.max_log_moneyness(),
+                derivative.time_to_maturity(),
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.gamma(
+                derivative.log_moneyness(),
+                None,
+                derivative.time_to_maturity(),
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.gamma(
+                derivative.log_moneyness(),
+                derivative.max_log_moneyness(),
+                None,
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.gamma(
+                derivative.log_moneyness(),
+                derivative.max_log_moneyness(),
+                derivative.time_to_maturity(),
+                None,
+            )
+
     def test_check_vega(self):
         m = BSAmericanBinaryOption()
 
@@ -108,6 +308,106 @@ def test_check_vega(self):
         expect = torch.tensor([0.0])
         assert_close(result, expect)
 
+    @pytest.mark.parametrize("call", [True])
+    def test_vega_2(self, call: bool):
+        m = BSAmericanBinaryOption(call=call)
+        with pytest.raises(ValueError):
+            m.vega(
+                torch.tensor(0.0),
+                torch.tensor(0.0),
+                torch.tensor(-1.0),
+                torch.tensor(0.2),
+            )
+        with pytest.raises(ValueError):
+            m.vega(
+                torch.tensor(0.0),
+                torch.tensor(0.0),
+                torch.tensor(1.0),
+                torch.tensor(-0.2),
+            )
+
+    @pytest.mark.parametrize("call", [True])
+    def test_vega_3(self, call: bool):
+        derivative = AmericanBinaryOption(BrownianStock(), call=call)
+        m = BSAmericanBinaryOption.from_derivative(derivative)
+        m2 = BSAmericanBinaryOption(call=call)
+        with pytest.raises(AttributeError):
+            m.vega(None, torch.tensor(1), torch.tensor(2), torch.tensor(3))
+        with pytest.raises(AttributeError):
+            m.vega(torch.tensor(1), None, torch.tensor(2), torch.tensor(3))
+        with pytest.raises(AttributeError):
+            m.vega(torch.tensor(1), torch.tensor(2), None, torch.tensor(3))
+        # ToDo: #530
+        # with pytest.raises(AttributeError):
+        #     m.vega(torch.tensor(1), torch.tensor(2), torch.tensor(3), None)
+        torch.manual_seed(42)
+        derivative.simulate(n_paths=1)
+        result = m.vega()
+        expect = m2.vega(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        # ToDo: [..., :-1] should be removed
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.vega(
+            None,
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.vega(
+            derivative.log_moneyness(),
+            None,
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.vega(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            None,
+            derivative.underlier.volatility,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.vega(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            None,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        with pytest.raises(ValueError):
+            m2.vega(
+                None,
+                derivative.max_log_moneyness(),
+                derivative.time_to_maturity(),
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.vega(
+                derivative.log_moneyness(),
+                None,
+                derivative.time_to_maturity(),
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.vega(
+                derivative.log_moneyness(),
+                derivative.max_log_moneyness(),
+                None,
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.vega(
+                derivative.log_moneyness(),
+                derivative.max_log_moneyness(),
+                derivative.time_to_maturity(),
+                None,
+            )
+
     def test_check_theta(self):
         m = BSAmericanBinaryOption()
 
@@ -131,6 +431,106 @@ def test_check_theta(self):
         expect = torch.tensor([0.0])
         assert_close(result, expect)
 
+    @pytest.mark.parametrize("call", [True])
+    def test_theta_2(self, call: bool):
+        m = BSAmericanBinaryOption(call=call)
+        with pytest.raises(ValueError):
+            m.theta(
+                torch.tensor(0.0),
+                torch.tensor(0.0),
+                torch.tensor(-1.0),
+                torch.tensor(0.2),
+            )
+        with pytest.raises(ValueError):
+            m.theta(
+                torch.tensor(0.0),
+                torch.tensor(0.0),
+                torch.tensor(1.0),
+                torch.tensor(-0.2),
+            )
+
+    @pytest.mark.parametrize("call", [True])
+    def test_theta_3(self, call: bool):
+        derivative = AmericanBinaryOption(BrownianStock(), call=call)
+        m = BSAmericanBinaryOption.from_derivative(derivative)
+        m2 = BSAmericanBinaryOption(call=call)
+        with pytest.raises(AttributeError):
+            m.theta(None, torch.tensor(1), torch.tensor(2), torch.tensor(3))
+        with pytest.raises(AttributeError):
+            m.theta(torch.tensor(1), None, torch.tensor(2), torch.tensor(3))
+        with pytest.raises(AttributeError):
+            m.theta(torch.tensor(1), torch.tensor(2), None, torch.tensor(3))
+        # ToDo: #530
+        # with pytest.raises(AttributeError):
+        #     m.theta(torch.tensor(1), torch.tensor(2), torch.tensor(3), None)
+        torch.manual_seed(42)
+        derivative.simulate(n_paths=1)
+        result = m.theta()
+        expect = m2.theta(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        # ToDo: [..., :-1] should be removed
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.theta(
+            None,
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.theta(
+            derivative.log_moneyness(),
+            None,
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.theta(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            None,
+            derivative.underlier.volatility,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.theta(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            None,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        with pytest.raises(ValueError):
+            m2.theta(
+                None,
+                derivative.max_log_moneyness(),
+                derivative.time_to_maturity(),
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.theta(
+                derivative.log_moneyness(),
+                None,
+                derivative.time_to_maturity(),
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.theta(
+                derivative.log_moneyness(),
+                derivative.max_log_moneyness(),
+                None,
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.theta(
+                derivative.log_moneyness(),
+                derivative.max_log_moneyness(),
+                derivative.time_to_maturity(),
+                None,
+            )
+
     def test_check_price(self):
         m = BSAmericanBinaryOption()
 
@@ -154,6 +554,105 @@ def test_check_price(self):
         expect = torch.tensor([0.0])
         assert_close(result, expect)
 
+    @pytest.mark.parametrize("call", [True])
+    def test_price_3(self, call: bool):
+        m = BSAmericanBinaryOption(call=call)
+        with pytest.raises(ValueError):
+            m.price(
+                torch.tensor(0.0),
+                torch.tensor(0.0),
+                torch.tensor(-1.0),
+                torch.tensor(0.2),
+            )
+        with pytest.raises(ValueError):
+            m.price(
+                torch.tensor(0.0),
+                torch.tensor(0.0),
+                torch.tensor(1.0),
+                torch.tensor(-0.2),
+            )
+
+    @pytest.mark.parametrize("call", [True])
+    def test_price_4(self, call: bool):
+        derivative = AmericanBinaryOption(BrownianStock(), call=call)
+        m = BSAmericanBinaryOption.from_derivative(derivative)
+        m2 = BSAmericanBinaryOption(call=call)
+        with pytest.raises(AttributeError):
+            m.price(None, torch.tensor(1), torch.tensor(2), torch.tensor(3))
+        with pytest.raises(AttributeError):
+            m.price(torch.tensor(1), None, torch.tensor(2), torch.tensor(3))
+        with pytest.raises(AttributeError):
+            m.price(torch.tensor(1), torch.tensor(2), None, torch.tensor(3))
+        # ToDo: #530
+        # with pytest.raises(AttributeError):
+        #     m.price(torch.tensor(1), torch.tensor(2), torch.tensor(3), None)
+        torch.manual_seed(42)
+        derivative.simulate(n_paths=1)
+        result = m.price()
+        expect = m2.price(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result, expect)
+        result = m.price(
+            None,
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result, expect)
+        result = m.price(
+            derivative.log_moneyness(),
+            None,
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result, expect)
+        result = m.price(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            None,
+            derivative.underlier.volatility,
+        )
+        assert_close(result, expect)
+        result = m.price(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            None,
+        )
+        assert_close(result, expect)
+        with pytest.raises(ValueError):
+            m2.price(
+                None,
+                derivative.max_log_moneyness(),
+                derivative.time_to_maturity(),
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.price(
+                derivative.log_moneyness(),
+                None,
+                derivative.time_to_maturity(),
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.price(
+                derivative.log_moneyness(),
+                derivative.max_log_moneyness(),
+                None,
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.price(
+                derivative.log_moneyness(),
+                derivative.max_log_moneyness(),
+                derivative.time_to_maturity(),
+                None,
+            )
+
     def test_check_price_monte_carlo(self):
         torch.manual_seed(42)
 
@@ -198,6 +697,27 @@ def test_vega_and_gamma(self):
         gamma = m.gamma(spot.log(), spot.log(), t, v)
         assert_close(vega, spot.square() * v * t * gamma, atol=1e-3, rtol=0)
 
+    @pytest.mark.parametrize("call", [True])
+    def test_vega_and_gamma_2(self, call: bool):
+        derivative = AmericanBinaryOption(BrownianStock(), call=call)
+        m = BSAmericanBinaryOption.from_derivative(derivative)
+        torch.manual_seed(42)
+        derivative.simulate(n_paths=1)
+        vega = m.vega()
+        gamma = m.gamma()
+        # ToDo: [..., :-1] should be removed
+        assert_close(
+            vega[..., :-1],
+            (
+                derivative.underlier.spot.square()
+                * derivative.underlier.volatility
+                * derivative.time_to_maturity()
+                * gamma
+            )[..., :-1],
+            atol=1e-3,
+            rtol=0,
+        )
+
     def test_features(self):
         m = BSAmericanBinaryOption()
         assert m.inputs() == [
@@ -226,6 +746,81 @@ def test_implied_volatility(self):
         expect = input[:, -1]
         assert_close(result, expect, atol=1e-4, rtol=1e-4, check_stride=False)
 
+    @pytest.mark.parametrize("call", [True])
+    def test_implied_volatility_2(self, call: bool):
+        derivative = AmericanBinaryOption(BrownianStock(), call=call)
+        m = BSAmericanBinaryOption.from_derivative(derivative)
+        m2 = BSAmericanBinaryOption(call=call)
+        with pytest.raises(AttributeError):
+            m.implied_volatility(
+                None, torch.tensor(1), torch.tensor(2), torch.tensor(3)
+            )
+        with pytest.raises(AttributeError):
+            m.implied_volatility(
+                torch.tensor(1), None, torch.tensor(2), torch.tensor(3)
+            )
+        with pytest.raises(AttributeError):
+            m.implied_volatility(
+                torch.tensor(1), torch.tensor(2), None, torch.tensor(3)
+            )
+        # ToDo: #530
+        # with pytest.raises(AttributeError):
+        #     m.implied_volatility(torch.tensor(1), torch.tensor(2), torch.tensor(3), None)
+        torch.manual_seed(42)
+        derivative.simulate(n_paths=1)
+        with pytest.raises(ValueError):
+            m.implied_volatility()
+        result = m.implied_volatility(price=derivative.underlier.spot)
+        expect = m2.implied_volatility(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.spot,
+        )
+        assert_close(result, expect)
+        result = m.implied_volatility(
+            None,
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.spot,
+        )
+        assert_close(result, expect)
+        result = m.implied_volatility(
+            derivative.log_moneyness(),
+            None,
+            derivative.time_to_maturity(),
+            derivative.underlier.spot,
+        )
+        assert_close(result, expect)
+        result = m.implied_volatility(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            None,
+            derivative.underlier.spot,
+        )
+        assert_close(result, expect)
+        with pytest.raises(ValueError):
+            m2.implied_volatility(
+                None,
+                derivative.max_log_moneyness(),
+                derivative.time_to_maturity(),
+                derivative.underlier.spot,
+            )
+        with pytest.raises(ValueError):
+            m2.implied_volatility(
+                derivative.log_moneyness(),
+                None,
+                derivative.time_to_maturity(),
+                derivative.underlier.spot,
+            )
+        with pytest.raises(ValueError):
+            m2.implied_volatility(
+                derivative.log_moneyness(),
+                derivative.max_log_moneyness(),
+                None,
+                derivative.underlier.spot,
+            )
+
     def test_shape(self):
         torch.distributions.Distribution.set_default_validate_args(False)
 
diff --git a/tests/nn/modules/bs/test_european.py b/tests/nn/modules/bs/test_european.py
index d93b01f4..5b4e253d 100644
--- a/tests/nn/modules/bs/test_european.py
+++ b/tests/nn/modules/bs/test_european.py
@@ -230,6 +230,48 @@ def test_delta_3(self, call: bool):
         expect = torch.full_like(result, 0.5 if call else -0.5)
         assert_close(result, expect)
 
+    @pytest.mark.parametrize("call", [True, False])
+    def test_delta_4(self, call: bool):
+        derivative = EuropeanOption(BrownianStock(), call=call)
+        m = BSEuropeanOption.from_derivative(derivative)
+        m2 = BSEuropeanOption(call=call)
+        with pytest.raises(AttributeError):
+            m.delta(None, torch.tensor(1), torch.tensor(2))
+        with pytest.raises(AttributeError):
+            m.delta(torch.tensor(1), None, torch.tensor(2))
+        # ToDo: #530
+        # with pytest.raises(AttributeError):
+        #     m.delta(torch.tensor(1), torch.tensor(2), None)
+        torch.manual_seed(42)
+        derivative.simulate(n_paths=1)
+        result = m.delta()
+        expect = m2.delta(
+            derivative.log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result, expect)
+        result = m.delta(
+            None, derivative.time_to_maturity(), derivative.underlier.volatility
+        )
+        assert_close(result, expect)
+        result = m.delta(
+            derivative.log_moneyness(), None, derivative.underlier.volatility
+        )
+        assert_close(result, expect)
+        result = m.delta(
+            derivative.log_moneyness(), derivative.time_to_maturity(), None
+        )
+        assert_close(result, expect)
+        with pytest.raises(ValueError):
+            m2.delta(
+                None, derivative.time_to_maturity(), derivative.underlier.volatility
+            )
+        with pytest.raises(ValueError):
+            m2.delta(derivative.log_moneyness(), None, derivative.underlier.volatility)
+        with pytest.raises(ValueError):
+            m2.delta(derivative.log_moneyness(), derivative.time_to_maturity(), None)
+
     def test_gamma_1(self):
         m = BSEuropeanOption()
         result = m.gamma(torch.tensor(0.0), torch.tensor(1.0), torch.tensor(0.2))
@@ -265,6 +307,48 @@ def test_gamma_3(self, call: bool):
         expect = torch.full_like(result, float("inf"))
         assert_close(result, expect)
 
+    @pytest.mark.parametrize("call", [True, False])
+    def test_gamma_4(self, call: bool):
+        derivative = EuropeanOption(BrownianStock(), call=call)
+        m = BSEuropeanOption.from_derivative(derivative)
+        m2 = BSEuropeanOption(call=call)
+        with pytest.raises(AttributeError):
+            m.gamma(None, torch.tensor(1), torch.tensor(2))
+        with pytest.raises(AttributeError):
+            m.gamma(torch.tensor(1), None, torch.tensor(2))
+        # ToDo: #530
+        # with pytest.raises(AttributeError):
+        #     m.gamma(torch.tensor(1), torch.tensor(2), None)
+        torch.manual_seed(42)
+        derivative.simulate(n_paths=1)
+        result = m.gamma()
+        expect = m2.gamma(
+            derivative.log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result, expect)
+        result = m.gamma(
+            None, derivative.time_to_maturity(), derivative.underlier.volatility
+        )
+        assert_close(result, expect)
+        result = m.gamma(
+            derivative.log_moneyness(), None, derivative.underlier.volatility
+        )
+        assert_close(result, expect)
+        result = m.gamma(
+            derivative.log_moneyness(), derivative.time_to_maturity(), None
+        )
+        assert_close(result, expect)
+        with pytest.raises(ValueError):
+            m2.gamma(
+                None, derivative.time_to_maturity(), derivative.underlier.volatility
+            )
+        with pytest.raises(ValueError):
+            m2.gamma(derivative.log_moneyness(), None, derivative.underlier.volatility)
+        with pytest.raises(ValueError):
+            m2.gamma(derivative.log_moneyness(), derivative.time_to_maturity(), None)
+
     def test_price_1(self):
         m = BSEuropeanOption()
         result = m.price(0.0, 1.0, 0.2)
@@ -306,6 +390,48 @@ def test_price_3(self, call: bool):
         expect = torch.full_like(result, 0.0)
         assert_close(result, expect)
 
+    @pytest.mark.parametrize("call", [True, False])
+    def test_price_4(self, call: bool):
+        derivative = EuropeanOption(BrownianStock(), call=call)
+        m = BSEuropeanOption.from_derivative(derivative)
+        m2 = BSEuropeanOption(call=call)
+        with pytest.raises(AttributeError):
+            m.price(None, torch.tensor(1), torch.tensor(2))
+        with pytest.raises(AttributeError):
+            m.price(torch.tensor(1), None, torch.tensor(2))
+        # ToDo: #530
+        # with pytest.raises(AttributeError):
+        #     m.price(torch.tensor(1), torch.tensor(2), None)
+        torch.manual_seed(42)
+        derivative.simulate(n_paths=1)
+        result = m.price()
+        expect = m2.price(
+            derivative.log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result, expect)
+        result = m.price(
+            None, derivative.time_to_maturity(), derivative.underlier.volatility
+        )
+        assert_close(result, expect)
+        result = m.price(
+            derivative.log_moneyness(), None, derivative.underlier.volatility
+        )
+        assert_close(result, expect)
+        result = m.price(
+            derivative.log_moneyness(), derivative.time_to_maturity(), None
+        )
+        assert_close(result, expect)
+        with pytest.raises(ValueError):
+            m2.price(
+                None, derivative.time_to_maturity(), derivative.underlier.volatility
+            )
+        with pytest.raises(ValueError):
+            m2.price(derivative.log_moneyness(), None, derivative.underlier.volatility)
+        with pytest.raises(ValueError):
+            m2.price(derivative.log_moneyness(), derivative.time_to_maturity(), None)
+
     def test_implied_volatility(self):
         input = torch.tensor([[0.0, 0.1, 0.01], [0.0, 0.1, 0.02], [0.0, 0.1, 0.03]])
         m = BSEuropeanOption()
@@ -315,6 +441,46 @@ def test_implied_volatility(self):
         expect = input[:, 2]
         assert_close(result, expect, check_stride=False)
 
+    @pytest.mark.parametrize("call", [True, False])
+    def test_implied_volatility_2(self, call: bool):
+        derivative = EuropeanOption(BrownianStock(), call=call)
+        m = BSEuropeanOption.from_derivative(derivative)
+        m2 = BSEuropeanOption(call=call)
+        with pytest.raises(AttributeError):
+            m.implied_volatility(None, torch.tensor(1), torch.tensor(2))
+        with pytest.raises(AttributeError):
+            m.implied_volatility(torch.tensor(1), None, torch.tensor(2))
+        # ToDo: #530
+        # with pytest.raises(AttributeError):
+        #     m.implied_volatility(torch.tensor(1), torch.tensor(2), None)
+        torch.manual_seed(42)
+        derivative.simulate(n_paths=1)
+        with pytest.raises(ValueError):
+            m.implied_volatility()
+        result = m.implied_volatility(price=derivative.underlier.spot)
+        expect = m2.implied_volatility(
+            derivative.log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.spot,
+        )
+        assert_close(result, expect)
+        result = m.implied_volatility(
+            None, derivative.time_to_maturity(), derivative.underlier.spot
+        )
+        assert_close(result, expect)
+        result = m.implied_volatility(
+            derivative.log_moneyness(), None, derivative.underlier.spot
+        )
+        assert_close(result, expect)
+        with pytest.raises(ValueError):
+            m2.implied_volatility(
+                None, derivative.time_to_maturity(), derivative.underlier.spot
+            )
+        with pytest.raises(ValueError):
+            m2.implied_volatility(
+                derivative.log_moneyness(), None, derivative.underlier.spot
+            )
+
     def test_vega(self):
         input = torch.tensor([[0.0, 0.1, 0.2], [0.0, 0.2, 0.2], [0.0, 0.3, 0.2]])
         m = BSEuropeanOption()
@@ -349,6 +515,46 @@ def test_vega_2(self, call: bool):
         expect = torch.full_like(result, 0.0)
         assert_close(result, expect)
 
+    @pytest.mark.parametrize("call", [True, False])
+    def test_vega_3(self, call: bool):
+        derivative = EuropeanOption(BrownianStock(), call=call)
+        m = BSEuropeanOption.from_derivative(derivative)
+        m2 = BSEuropeanOption(call=call)
+        with pytest.raises(AttributeError):
+            m.vega(None, torch.tensor(1), torch.tensor(2))
+        with pytest.raises(AttributeError):
+            m.vega(torch.tensor(1), None, torch.tensor(2))
+        # ToDo: #530
+        # with pytest.raises(AttributeError):
+        #     m.vega(torch.tensor(1), torch.tensor(2), None)
+        torch.manual_seed(42)
+        derivative.simulate(n_paths=1)
+        result = m.vega()
+        expect = m2.vega(
+            derivative.log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result, expect)
+        result = m.vega(
+            None, derivative.time_to_maturity(), derivative.underlier.volatility
+        )
+        assert_close(result, expect)
+        result = m.vega(
+            derivative.log_moneyness(), None, derivative.underlier.volatility
+        )
+        assert_close(result, expect)
+        result = m.vega(derivative.log_moneyness(), derivative.time_to_maturity(), None)
+        assert_close(result, expect)
+        with pytest.raises(ValueError):
+            m2.vega(
+                None, derivative.time_to_maturity(), derivative.underlier.volatility
+            )
+        with pytest.raises(ValueError):
+            m2.vega(derivative.log_moneyness(), None, derivative.underlier.volatility)
+        with pytest.raises(ValueError):
+            m2.vega(derivative.log_moneyness(), derivative.time_to_maturity(), None)
+
     def test_vega_and_gamma(self):
         m = BSEuropeanOption()
         # vega = spot^2 * sigma * (T - t) * gamma
@@ -360,6 +566,24 @@ def test_vega_and_gamma(self):
         gamma = m.gamma(spot.log(), t, v)
         assert_close(vega, spot.square() * v * t * gamma, atol=1e-3, rtol=0)
 
+    @pytest.mark.parametrize("call", [True, False])
+    def test_vega_and_gamma_2(self, call: bool):
+        derivative = EuropeanOption(BrownianStock(), call=call)
+        m = BSEuropeanOption.from_derivative(derivative)
+        torch.manual_seed(42)
+        derivative.simulate(n_paths=1)
+        vega = m.vega()
+        gamma = m.gamma()
+        assert_close(
+            vega,
+            derivative.underlier.spot.square()
+            * derivative.underlier.volatility
+            * derivative.time_to_maturity()
+            * gamma,
+            atol=1e-3,
+            rtol=0,
+        )
+
     def test_theta(self):
         input = torch.tensor([[0.0, 0.1, 0.2], [0.0, 0.2, 0.2], [0.0, 0.3, 0.2]])
         m = BSEuropeanOption(strike=100)
@@ -394,6 +618,48 @@ def test_theta_2(self, call: bool):
         expect = torch.full_like(result, -0.0)
         assert_close(result, expect)
 
+    @pytest.mark.parametrize("call", [True, False])
+    def test_theta_3(self, call: bool):
+        derivative = EuropeanOption(BrownianStock(), call=call)
+        m = BSEuropeanOption.from_derivative(derivative)
+        m2 = BSEuropeanOption(call=call)
+        with pytest.raises(AttributeError):
+            m.theta(None, torch.tensor(1), torch.tensor(2))
+        with pytest.raises(AttributeError):
+            m.theta(torch.tensor(1), None, torch.tensor(2))
+        # ToDo: #530
+        # with pytest.raises(AttributeError):
+        #     m.theta(torch.tensor(1), torch.tensor(2), None)
+        torch.manual_seed(42)
+        derivative.simulate(n_paths=1)
+        result = m.theta()
+        expect = m2.theta(
+            derivative.log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result, expect)
+        result = m.theta(
+            None, derivative.time_to_maturity(), derivative.underlier.volatility
+        )
+        assert_close(result, expect)
+        result = m.theta(
+            derivative.log_moneyness(), None, derivative.underlier.volatility
+        )
+        assert_close(result, expect)
+        result = m.theta(
+            derivative.log_moneyness(), derivative.time_to_maturity(), None
+        )
+        assert_close(result, expect)
+        with pytest.raises(ValueError):
+            m2.theta(
+                None, derivative.time_to_maturity(), derivative.underlier.volatility
+            )
+        with pytest.raises(ValueError):
+            m2.theta(derivative.log_moneyness(), None, derivative.underlier.volatility)
+        with pytest.raises(ValueError):
+            m2.theta(derivative.log_moneyness(), derivative.time_to_maturity(), None)
+
     def test_example(self):
         from pfhedge.instruments import BrownianStock
         from pfhedge.instruments import EuropeanOption
diff --git a/tests/nn/modules/bs/test_european_binary.py b/tests/nn/modules/bs/test_european_binary.py
index f822195d..90d391d5 100644
--- a/tests/nn/modules/bs/test_european_binary.py
+++ b/tests/nn/modules/bs/test_european_binary.py
@@ -1,6 +1,5 @@
 import pytest
 import torch
-from torch import Tensor
 from torch.testing import assert_close
 
 from pfhedge.features._getter import get_feature
@@ -185,18 +184,183 @@ def test_delta(self):
         expect = torch.tensor(6.3047)
         assert_close(result, expect, atol=1e-4, rtol=1e-4)
 
+    @pytest.mark.parametrize("call", [True, False])
+    def test_delta_2(self, call: bool):
+        m = BSEuropeanBinaryOption(call=call)
+        with pytest.raises(ValueError):
+            m.delta(torch.tensor(0.0), torch.tensor(-1.0), torch.tensor(0.2))
+        with pytest.raises(ValueError):
+            m.delta(torch.tensor(0.0), torch.tensor(1.0), torch.tensor(-0.2))
+        result = m.delta(torch.tensor(1.0), torch.tensor(1.0), torch.tensor(0))
+        expect = torch.full_like(result, 0.0)
+        assert_close(result, expect)
+        result = m.delta(torch.tensor(1.0), torch.tensor(0.0), torch.tensor(0.1))
+        expect = torch.full_like(result, 0.0)
+        assert_close(result, expect)
+        result = m.delta(torch.tensor(1.0), torch.tensor(0.0), torch.tensor(0.0))
+        expect = torch.full_like(result, 0.0)
+        assert_close(result, expect)
+        result = m.delta(torch.tensor(0.0), torch.tensor(0.0), torch.tensor(0.1))
+        expect = torch.full_like(result, float("inf") if call else -float("inf"))
+        assert_close(result, expect)
+        result = m.delta(torch.tensor(0.0), torch.tensor(0.0), torch.tensor(0.0))
+        expect = torch.full_like(result, float("inf") if call else -float("inf"))
+        assert_close(result, expect)
+
+    @pytest.mark.parametrize("call", [True, False])
+    def test_delta_3(self, call: bool):
+        derivative = EuropeanBinaryOption(BrownianStock(), call=call)
+        m = BSEuropeanBinaryOption.from_derivative(derivative)
+        m2 = BSEuropeanBinaryOption(call=call)
+        with pytest.raises(AttributeError):
+            m.delta(None, torch.tensor(1), torch.tensor(2))
+        with pytest.raises(AttributeError):
+            m.delta(torch.tensor(1), None, torch.tensor(2))
+        # ToDo: #530
+        # with pytest.raises(AttributeError):
+        #     m.delta(torch.tensor(1), torch.tensor(2), None)
+        torch.manual_seed(42)
+        derivative.simulate(n_paths=1)
+        result = m.delta()
+        expect = m2.delta(
+            derivative.log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result, expect)
+        result = m.delta(
+            None, derivative.time_to_maturity(), derivative.underlier.volatility
+        )
+        assert_close(result, expect)
+        result = m.delta(
+            derivative.log_moneyness(), None, derivative.underlier.volatility
+        )
+        assert_close(result, expect)
+        result = m.delta(
+            derivative.log_moneyness(), derivative.time_to_maturity(), None
+        )
+        assert_close(result, expect)
+        with pytest.raises(ValueError):
+            m2.delta(
+                None, derivative.time_to_maturity(), derivative.underlier.volatility
+            )
+        with pytest.raises(ValueError):
+            m2.delta(derivative.log_moneyness(), None, derivative.underlier.volatility)
+        with pytest.raises(ValueError):
+            m2.delta(derivative.log_moneyness(), derivative.time_to_maturity(), None)
+
     def test_gamma(self):
         m = BSEuropeanBinaryOption()
         result = m.gamma(torch.tensor(0.01), torch.tensor(1.0), torch.tensor(0.2))
         expect = torch.tensor(-1.4645787477493286)
         assert_close(result, expect)
 
+    @pytest.mark.parametrize("call", [True, False])
+    def test_gamma_2(self, call: bool):
+        m = BSEuropeanBinaryOption(call=call)
+        with pytest.raises(ValueError):
+            m.gamma(torch.tensor(0.0), torch.tensor(-1.0), torch.tensor(0.2))
+        with pytest.raises(ValueError):
+            m.gamma(torch.tensor(0.0), torch.tensor(1.0), torch.tensor(-0.2))
+
+    @pytest.mark.parametrize("call", [True, False])
+    def test_gamma_3(self, call: bool):
+        derivative = EuropeanBinaryOption(BrownianStock(), call=call)
+        m = BSEuropeanBinaryOption.from_derivative(derivative)
+        m2 = BSEuropeanBinaryOption(call=call)
+        with pytest.raises(AttributeError):
+            m.gamma(None, torch.tensor(1), torch.tensor(2))
+        with pytest.raises(AttributeError):
+            m.gamma(torch.tensor(1), None, torch.tensor(2))
+        # ToDo: #530
+        # with pytest.raises(AttributeError):
+        #     m.gamma(torch.tensor(1), torch.tensor(2), None)
+        torch.manual_seed(42)
+        derivative.simulate(n_paths=1)
+        result = m.gamma()
+        expect = m2.gamma(
+            derivative.log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        # ToDo: [..., :-1] should be removed
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.gamma(
+            None, derivative.time_to_maturity(), derivative.underlier.volatility
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.gamma(
+            derivative.log_moneyness(), None, derivative.underlier.volatility
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.gamma(
+            derivative.log_moneyness(), derivative.time_to_maturity(), None
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        with pytest.raises(ValueError):
+            m2.gamma(
+                None, derivative.time_to_maturity(), derivative.underlier.volatility
+            )
+        with pytest.raises(ValueError):
+            m2.gamma(derivative.log_moneyness(), None, derivative.underlier.volatility)
+        with pytest.raises(ValueError):
+            m2.gamma(derivative.log_moneyness(), derivative.time_to_maturity(), None)
+
     def test_vega(self):
         m = BSEuropeanBinaryOption()
         result = m.vega(torch.tensor(0.0), torch.tensor(0.1), torch.tensor(0.2))
         expect = torch.tensor(-0.06305)
         assert_close(result, expect, atol=1e-4, rtol=1e-4)
 
+    @pytest.mark.parametrize("call", [True, False])
+    def test_vega_2(self, call: bool):
+        m = BSEuropeanBinaryOption(call=call)
+        with pytest.raises(ValueError):
+            m.vega(torch.tensor(0.0), torch.tensor(-1.0), torch.tensor(0.2))
+        with pytest.raises(ValueError):
+            m.vega(torch.tensor(0.0), torch.tensor(1.0), torch.tensor(-0.2))
+
+    @pytest.mark.parametrize("call", [True, False])
+    def test_vega_3(self, call: bool):
+        derivative = EuropeanBinaryOption(BrownianStock(), call=call)
+        m = BSEuropeanBinaryOption.from_derivative(derivative)
+        m2 = BSEuropeanBinaryOption(call=call)
+        with pytest.raises(AttributeError):
+            m.vega(None, torch.tensor(1), torch.tensor(2))
+        with pytest.raises(AttributeError):
+            m.vega(torch.tensor(1), None, torch.tensor(2))
+        # ToDo: #530
+        # with pytest.raises(AttributeError):
+        #     m.vega(torch.tensor(1), torch.tensor(2), None)
+        torch.manual_seed(42)
+        derivative.simulate(n_paths=1)
+        result = m.vega()
+        expect = m2.vega(
+            derivative.log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        # ToDo: [..., :-1] should be removed
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.vega(
+            None, derivative.time_to_maturity(), derivative.underlier.volatility
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.vega(
+            derivative.log_moneyness(), None, derivative.underlier.volatility
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.vega(derivative.log_moneyness(), derivative.time_to_maturity(), None)
+        assert_close(result[..., :-1], expect[..., :-1])
+        with pytest.raises(ValueError):
+            m2.vega(
+                None, derivative.time_to_maturity(), derivative.underlier.volatility
+            )
+        with pytest.raises(ValueError):
+            m2.vega(derivative.log_moneyness(), None, derivative.underlier.volatility)
+        with pytest.raises(ValueError):
+            m2.vega(derivative.log_moneyness(), derivative.time_to_maturity(), None)
+
     def test_vega_and_gamma(self):
         m = BSEuropeanBinaryOption()
         # vega = spot^2 * sigma * (T - t) * gamma
@@ -208,12 +372,84 @@ def test_vega_and_gamma(self):
         gamma = m.gamma(spot.log(), t, v)
         assert_close(vega, spot.square() * v * t * gamma, atol=1e-3, rtol=0)
 
+    @pytest.mark.parametrize("call", [True, False])
+    def test_vega_and_gamma_2(self, call: bool):
+        derivative = EuropeanBinaryOption(BrownianStock(), call=call)
+        m = BSEuropeanBinaryOption.from_derivative(derivative)
+        torch.manual_seed(42)
+        derivative.simulate(n_paths=1)
+        vega = m.vega()
+        gamma = m.gamma()
+        # ToDo: [..., :-1] should be removed
+        assert_close(
+            vega[..., :-1],
+            (
+                derivative.underlier.spot.square()
+                * derivative.underlier.volatility
+                * derivative.time_to_maturity()
+                * gamma
+            )[..., :-1],
+            atol=1e-3,
+            rtol=0,
+        )
+
     def test_theta(self):
         m = BSEuropeanBinaryOption()
         result = m.theta(torch.tensor(0.0), torch.tensor(0.1), torch.tensor(0.2))
         expect = torch.tensor(0.0630)
         assert_close(result, expect, atol=1e-4, rtol=1e-4)
 
+    @pytest.mark.parametrize("call", [True, False])
+    def test_theta_2(self, call: bool):
+        m = BSEuropeanBinaryOption(call=call)
+        with pytest.raises(ValueError):
+            m.theta(torch.tensor(0.0), torch.tensor(-1.0), torch.tensor(0.2))
+        with pytest.raises(ValueError):
+            m.theta(torch.tensor(0.0), torch.tensor(1.0), torch.tensor(-0.2))
+
+    @pytest.mark.parametrize("call", [True, False])
+    def test_theta_3(self, call: bool):
+        derivative = EuropeanBinaryOption(BrownianStock(), call=call)
+        m = BSEuropeanBinaryOption.from_derivative(derivative)
+        m2 = BSEuropeanBinaryOption(call=call)
+        with pytest.raises(AttributeError):
+            m.theta(None, torch.tensor(1), torch.tensor(2))
+        with pytest.raises(AttributeError):
+            m.theta(torch.tensor(1), None, torch.tensor(2))
+        # ToDo: #530
+        # with pytest.raises(AttributeError):
+        #     m.theta(torch.tensor(1), torch.tensor(2), None)
+        torch.manual_seed(42)
+        derivative.simulate(n_paths=1)
+        result = m.theta()
+        expect = m2.theta(
+            derivative.log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        # ToDo: [..., :-1] should be removed
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.theta(
+            None, derivative.time_to_maturity(), derivative.underlier.volatility
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.theta(
+            derivative.log_moneyness(), None, derivative.underlier.volatility
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.theta(
+            derivative.log_moneyness(), derivative.time_to_maturity(), None
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        with pytest.raises(ValueError):
+            m2.theta(
+                None, derivative.time_to_maturity(), derivative.underlier.volatility
+            )
+        with pytest.raises(ValueError):
+            m2.theta(derivative.log_moneyness(), None, derivative.underlier.volatility)
+        with pytest.raises(ValueError):
+            m2.theta(derivative.log_moneyness(), derivative.time_to_maturity(), None)
+
     def test_price(self):
         m = BSEuropeanBinaryOption()
 
@@ -225,6 +461,77 @@ def test_price(self):
         expect = torch.tensor(0.4880)
         assert_close(result, expect, atol=1e-4, rtol=1e-4)
 
+    @pytest.mark.parametrize("call", [True, False])
+    def test_price_2(self, call: bool):
+        m = BSEuropeanBinaryOption(call=call)
+        with pytest.raises(ValueError):
+            m.price(torch.tensor(0.0), torch.tensor(-1.0), torch.tensor(0.2))
+        with pytest.raises(ValueError):
+            m.price(torch.tensor(0.0), torch.tensor(1.0), torch.tensor(-0.2))
+        result = m.price(
+            torch.tensor(1.0 if call else -1.0), torch.tensor(1.0), torch.tensor(0)
+        )
+        expect = torch.full_like(result, 1)
+        assert_close(result, expect)
+        result = m.price(
+            torch.tensor(1.0 if call else -1.0), torch.tensor(0.0), torch.tensor(0.1)
+        )
+        expect = torch.full_like(result, 1)
+        assert_close(result, expect)
+        result = m.price(
+            torch.tensor(1.0 if call else -1.0), torch.tensor(0.0), torch.tensor(0.0)
+        )
+        expect = torch.full_like(result, 1)
+        assert_close(result, expect)
+        result = m.price(torch.tensor(0.0), torch.tensor(0.0), torch.tensor(0.1))
+        expect = torch.full_like(result, 0.5)
+        assert_close(result, expect)
+        result = m.price(torch.tensor(0.0), torch.tensor(0.0), torch.tensor(0.0))
+        expect = torch.full_like(result, 0.5)
+        assert_close(result, expect)
+
+    @pytest.mark.parametrize("call", [True, False])
+    def test_price_3(self, call: bool):
+        derivative = EuropeanBinaryOption(BrownianStock(), call=call)
+        m = BSEuropeanBinaryOption.from_derivative(derivative)
+        m2 = BSEuropeanBinaryOption(call=call)
+        with pytest.raises(AttributeError):
+            m.price(None, torch.tensor(1), torch.tensor(2))
+        with pytest.raises(AttributeError):
+            m.price(torch.tensor(1), None, torch.tensor(2))
+        # ToDo: #530
+        # with pytest.raises(AttributeError):
+        #     m.price(torch.tensor(1), torch.tensor(2), None)
+        torch.manual_seed(42)
+        derivative.simulate(n_paths=1)
+        result = m.price()
+        expect = m2.price(
+            derivative.log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result, expect)
+        result = m.price(
+            None, derivative.time_to_maturity(), derivative.underlier.volatility
+        )
+        assert_close(result, expect)
+        result = m.price(
+            derivative.log_moneyness(), None, derivative.underlier.volatility
+        )
+        assert_close(result, expect)
+        result = m.price(
+            derivative.log_moneyness(), derivative.time_to_maturity(), None
+        )
+        assert_close(result, expect)
+        with pytest.raises(ValueError):
+            m2.price(
+                None, derivative.time_to_maturity(), derivative.underlier.volatility
+            )
+        with pytest.raises(ValueError):
+            m2.price(derivative.log_moneyness(), None, derivative.underlier.volatility)
+        with pytest.raises(ValueError):
+            m2.price(derivative.log_moneyness(), derivative.time_to_maturity(), None)
+
     def test_implied_volatility(self):
         lm = torch.full((3,), -0.01)
         t = torch.full((3,), 0.1)
@@ -236,6 +543,48 @@ def test_implied_volatility(self):
         result = BSEuropeanBinaryOption().price(lm, t, iv)
         assert_close(result, price, check_stride=False)
 
+    @pytest.mark.parametrize("call", [True, False])
+    def test_implied_volatility_2(self, call: bool):
+        derivative = EuropeanBinaryOption(BrownianStock(), call=call)
+        m = BSEuropeanBinaryOption.from_derivative(derivative)
+        m2 = BSEuropeanBinaryOption(call=call)
+        with pytest.raises(AttributeError):
+            m.implied_volatility()
+        with pytest.raises(AttributeError):
+            m.implied_volatility(None, torch.tensor(1), torch.tensor(1))
+        with pytest.raises(AttributeError):
+            m.implied_volatility(torch.tensor(1), None, torch.tensor(1))
+        # ToDo: #530
+        # with pytest.raises(AttributeError):
+        #     m.implied_volatility(torch.tensor(0), torch.tensor(0), None)
+        torch.manual_seed(42)
+        derivative.simulate(n_paths=1)
+        with pytest.raises(ValueError):
+            m.implied_volatility()
+        result = m.implied_volatility(price=derivative.underlier.spot)
+        expect = m2.implied_volatility(
+            derivative.log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.spot,
+        )
+        assert_close(result, expect)
+        result = m.implied_volatility(
+            None, derivative.time_to_maturity(), derivative.underlier.spot
+        )
+        assert_close(result, expect)
+        result = m.implied_volatility(
+            derivative.log_moneyness(), None, derivative.underlier.spot
+        )
+        assert_close(result, expect)
+        with pytest.raises(ValueError):
+            m2.implied_volatility(
+                None, derivative.time_to_maturity(), derivative.underlier.spot
+            )
+        with pytest.raises(ValueError):
+            m2.implied_volatility(
+                derivative.log_moneyness(), None, derivative.underlier.spot
+            )
+
     def test_example(self):
         from pfhedge.instruments import BrownianStock
         from pfhedge.instruments import EuropeanBinaryOption
diff --git a/tests/nn/modules/bs/test_lookback.py b/tests/nn/modules/bs/test_lookback.py
index 61690824..a66f20ca 100644
--- a/tests/nn/modules/bs/test_lookback.py
+++ b/tests/nn/modules/bs/test_lookback.py
@@ -66,6 +66,106 @@ def test_check_delta(self):
         expect = torch.tensor([0.0])
         assert_close(result, expect)
 
+    @pytest.mark.parametrize("call", [True])
+    def test_delta_2(self, call: bool):
+        m = BSLookbackOption(call=call)
+        with pytest.raises(ValueError):
+            m.delta(
+                torch.tensor(0.0),
+                torch.tensor(0.0),
+                torch.tensor(-1.0),
+                torch.tensor(0.2),
+            )
+        with pytest.raises(ValueError):
+            m.delta(
+                torch.tensor(0.0),
+                torch.tensor(0.0),
+                torch.tensor(1.0),
+                torch.tensor(-0.2),
+            )
+
+    @pytest.mark.parametrize("call", [True])
+    def test_delta_3(self, call: bool):
+        derivative = LookbackOption(BrownianStock(), call=call)
+        m = BSLookbackOption.from_derivative(derivative)
+        m2 = BSLookbackOption(call=call)
+        with pytest.raises(AttributeError):
+            m.delta(None, torch.zeros(1), torch.zeros(2), torch.zeros(3))
+        with pytest.raises(AttributeError):
+            m.delta(torch.zeros(1), None, torch.zeros(2), torch.zeros(3))
+        with pytest.raises(AttributeError):
+            m.delta(torch.zeros(1), torch.zeros(2), None, torch.zeros(3))
+        # ToDo: #530
+        # with pytest.raises(AttributeError):
+        #     m.delta(torch.zeros(1), torch.zeros(2), torch.zeros(3), None)
+        torch.manual_seed(42)
+        derivative.simulate(n_paths=1)
+        result = m.delta()
+        expect = m2.delta(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        # ToDo: [..., :-1] should be removed
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.delta(
+            None,
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.delta(
+            derivative.log_moneyness(),
+            None,
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.delta(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            None,
+            derivative.underlier.volatility,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.delta(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            None,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        with pytest.raises(ValueError):
+            m2.delta(
+                None,
+                derivative.max_log_moneyness(),
+                derivative.time_to_maturity(),
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.delta(
+                derivative.log_moneyness(),
+                None,
+                derivative.time_to_maturity(),
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.delta(
+                derivative.log_moneyness(),
+                derivative.max_log_moneyness(),
+                None,
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.delta(
+                derivative.log_moneyness(),
+                derivative.max_log_moneyness(),
+                derivative.time_to_maturity(),
+                None,
+            )
+
     def test_check_gamma(self):
         m = BSLookbackOption()
 
@@ -97,6 +197,106 @@ def test_check_gamma(self):
         expect = torch.tensor([0.0])
         assert_close(result, expect)
 
+    @pytest.mark.parametrize("call", [True])
+    def test_gamma_2(self, call: bool):
+        m = BSLookbackOption(call=call)
+        with pytest.raises(ValueError):
+            m.gamma(
+                torch.tensor(0.0),
+                torch.tensor(0.0),
+                torch.tensor(-1.0),
+                torch.tensor(0.2),
+            )
+        with pytest.raises(ValueError):
+            m.gamma(
+                torch.tensor(0.0),
+                torch.tensor(0.0),
+                torch.tensor(1.0),
+                torch.tensor(-0.2),
+            )
+
+    @pytest.mark.parametrize("call", [True])
+    def test_gamma_3(self, call: bool):
+        derivative = LookbackOption(BrownianStock(), call=call)
+        m = BSLookbackOption.from_derivative(derivative)
+        m2 = BSLookbackOption(call=call)
+        with pytest.raises(AttributeError):
+            m.gamma(None, torch.zeros(1), torch.zeros(2), torch.zeros(3))
+        with pytest.raises(AttributeError):
+            m.gamma(torch.zeros(1), None, torch.zeros(2), torch.zeros(3))
+        with pytest.raises(AttributeError):
+            m.gamma(torch.zeros(1), torch.zeros(2), None, torch.zeros(3))
+        # ToDo: #530
+        # with pytest.raises(AttributeError):
+        #     m.gamma(torch.zeros(1), torch.zeros(2), torch.zeros(3), None)
+        torch.manual_seed(42)
+        derivative.simulate(n_paths=1)
+        result = m.gamma()
+        expect = m2.gamma(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        # ToDo: [..., :-1] should be removed
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.gamma(
+            None,
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.gamma(
+            derivative.log_moneyness(),
+            None,
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.gamma(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            None,
+            derivative.underlier.volatility,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.gamma(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            None,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        with pytest.raises(ValueError):
+            m2.gamma(
+                None,
+                derivative.max_log_moneyness(),
+                derivative.time_to_maturity(),
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.gamma(
+                derivative.log_moneyness(),
+                None,
+                derivative.time_to_maturity(),
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.gamma(
+                derivative.log_moneyness(),
+                derivative.max_log_moneyness(),
+                None,
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.gamma(
+                derivative.log_moneyness(),
+                derivative.max_log_moneyness(),
+                derivative.time_to_maturity(),
+                None,
+            )
+
     def test_check_vega(self):
         m = BSLookbackOption()
 
@@ -128,6 +328,106 @@ def test_check_vega(self):
         expect = torch.tensor([0.0])
         assert_close(result, expect)
 
+    @pytest.mark.parametrize("call", [True])
+    def test_vega_2(self, call: bool):
+        m = BSLookbackOption(call=call)
+        with pytest.raises(ValueError):
+            m.vega(
+                torch.tensor(0.0),
+                torch.tensor(0.0),
+                torch.tensor(-1.0),
+                torch.tensor(0.2),
+            )
+        with pytest.raises(ValueError):
+            m.vega(
+                torch.tensor(0.0),
+                torch.tensor(0.0),
+                torch.tensor(1.0),
+                torch.tensor(-0.2),
+            )
+
+    @pytest.mark.parametrize("call", [True])
+    def test_vega_3(self, call: bool):
+        derivative = LookbackOption(BrownianStock(), call=call)
+        m = BSLookbackOption.from_derivative(derivative)
+        m2 = BSLookbackOption(call=call)
+        with pytest.raises(AttributeError):
+            m.vega(None, torch.tensor(1), torch.tensor(2), torch.tensor(3))
+        with pytest.raises(AttributeError):
+            m.vega(torch.tensor(1), None, torch.tensor(2), torch.tensor(3))
+        with pytest.raises(AttributeError):
+            m.vega(torch.tensor(1), torch.tensor(2), None, torch.tensor(3))
+        # ToDo: #530
+        # with pytest.raises(AttributeError):
+        #     m.vega(torch.tensor(1), torch.tensor(2), torch.tensor(3), None)
+        torch.manual_seed(42)
+        derivative.simulate(n_paths=1)
+        result = m.vega()
+        expect = m2.vega(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        # ToDo: [..., :-1] should be removed
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.vega(
+            None,
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.vega(
+            derivative.log_moneyness(),
+            None,
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.vega(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            None,
+            derivative.underlier.volatility,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.vega(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            None,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        with pytest.raises(ValueError):
+            m2.vega(
+                None,
+                derivative.max_log_moneyness(),
+                derivative.time_to_maturity(),
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.vega(
+                derivative.log_moneyness(),
+                None,
+                derivative.time_to_maturity(),
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.vega(
+                derivative.log_moneyness(),
+                derivative.max_log_moneyness(),
+                None,
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.vega(
+                derivative.log_moneyness(),
+                derivative.max_log_moneyness(),
+                derivative.time_to_maturity(),
+                None,
+            )
+
     def test_check_theta(self):
         m = BSLookbackOption()
 
@@ -159,6 +459,106 @@ def test_check_theta(self):
         expect = torch.tensor([0.0])
         assert_close(result, expect)
 
+    @pytest.mark.parametrize("call", [True])
+    def test_theta_2(self, call: bool):
+        m = BSLookbackOption(call=call)
+        with pytest.raises(ValueError):
+            m.theta(
+                torch.tensor(0.0),
+                torch.tensor(0.0),
+                torch.tensor(-1.0),
+                torch.tensor(0.2),
+            )
+        with pytest.raises(ValueError):
+            m.theta(
+                torch.tensor(0.0),
+                torch.tensor(0.0),
+                torch.tensor(1.0),
+                torch.tensor(-0.2),
+            )
+
+    @pytest.mark.parametrize("call", [True])
+    def test_theta_3(self, call: bool):
+        derivative = LookbackOption(BrownianStock(), call=call)
+        m = BSLookbackOption.from_derivative(derivative)
+        m2 = BSLookbackOption(call=call)
+        with pytest.raises(AttributeError):
+            m.theta(None, torch.tensor(1), torch.tensor(2), torch.tensor(3))
+        with pytest.raises(AttributeError):
+            m.theta(torch.tensor(1), None, torch.tensor(2), torch.tensor(3))
+        with pytest.raises(AttributeError):
+            m.theta(torch.tensor(1), torch.tensor(2), None, torch.tensor(3))
+        # ToDo: #530
+        # with pytest.raises(AttributeError):
+        #     m.theta(torch.tensor(1), torch.tensor(2), torch.tensor(3), None)
+        torch.manual_seed(42)
+        derivative.simulate(n_paths=1)
+        result = m.theta()
+        expect = m2.theta(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        # ToDo: [..., :-1] should be removed
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.theta(
+            None,
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.theta(
+            derivative.log_moneyness(),
+            None,
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.theta(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            None,
+            derivative.underlier.volatility,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.theta(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            None,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        with pytest.raises(ValueError):
+            m2.theta(
+                None,
+                derivative.max_log_moneyness(),
+                derivative.time_to_maturity(),
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.theta(
+                derivative.log_moneyness(),
+                None,
+                derivative.time_to_maturity(),
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.theta(
+                derivative.log_moneyness(),
+                derivative.max_log_moneyness(),
+                None,
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.theta(
+                derivative.log_moneyness(),
+                derivative.max_log_moneyness(),
+                derivative.time_to_maturity(),
+                None,
+            )
+
     def test_check_price(self):
         m = BSLookbackOption()
 
@@ -197,6 +597,106 @@ def test_check_price(self):
         result1 = compute_price(m, torch.tensor([[1e-5, 1e-5, 0.1, 0.2]]))
         assert_close(result0, result1, atol=1e-4, rtol=0)
 
+    @pytest.mark.parametrize("call", [True])
+    def test_price_3(self, call: bool):
+        m = BSLookbackOption(call=call)
+        with pytest.raises(ValueError):
+            m.price(
+                torch.tensor(0.0),
+                torch.tensor(0.0),
+                torch.tensor(-1.0),
+                torch.tensor(0.2),
+            )
+        with pytest.raises(ValueError):
+            m.price(
+                torch.tensor(0.0),
+                torch.tensor(0.0),
+                torch.tensor(1.0),
+                torch.tensor(-0.2),
+            )
+
+    @pytest.mark.parametrize("call", [True])
+    def test_price_4(self, call: bool):
+        derivative = LookbackOption(BrownianStock(), call=call)
+        m = BSLookbackOption.from_derivative(derivative)
+        m2 = BSLookbackOption(call=call)
+        with pytest.raises(AttributeError):
+            m.price(None, torch.tensor(1), torch.tensor(2), torch.tensor(3))
+        with pytest.raises(AttributeError):
+            m.price(torch.tensor(1), None, torch.tensor(2), torch.tensor(3))
+        with pytest.raises(AttributeError):
+            m.price(torch.tensor(1), torch.tensor(2), None, torch.tensor(3))
+        # ToDo: #530
+        # with pytest.raises(AttributeError):
+        #     m.price(torch.tensor(1), torch.tensor(2), torch.tensor(3), None)
+        torch.manual_seed(42)
+        derivative.simulate(n_paths=1)
+        result = m.price()
+        expect = m2.price(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        # ToDo: [..., :-1] should be removed
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.price(
+            None,
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.price(
+            derivative.log_moneyness(),
+            None,
+            derivative.time_to_maturity(),
+            derivative.underlier.volatility,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.price(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            None,
+            derivative.underlier.volatility,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        result = m.price(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            None,
+        )
+        assert_close(result[..., :-1], expect[..., :-1])
+        with pytest.raises(ValueError):
+            m2.price(
+                None,
+                derivative.max_log_moneyness(),
+                derivative.time_to_maturity(),
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.price(
+                derivative.log_moneyness(),
+                None,
+                derivative.time_to_maturity(),
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.price(
+                derivative.log_moneyness(),
+                derivative.max_log_moneyness(),
+                None,
+                derivative.underlier.volatility,
+            )
+        with pytest.raises(ValueError):
+            m2.price(
+                derivative.log_moneyness(),
+                derivative.max_log_moneyness(),
+                derivative.time_to_maturity(),
+                None,
+            )
+
     def test_check_price_monte_carlo(self):
         torch.manual_seed(42)
 
@@ -251,6 +751,68 @@ def test_implied_volatility(self):
         expect = input[:, -1]
         assert_close(result, expect, atol=1e-4, rtol=1e-4, check_stride=False)
 
+    @pytest.mark.parametrize("call", [True])
+    def test_implied_volatility_2(self, call: bool):
+        derivative = LookbackOption(BrownianStock(), call=call)
+        m = BSLookbackOption.from_derivative(derivative)
+        m2 = BSLookbackOption(call=call)
+        with pytest.raises(AttributeError):
+            m.implied_volatility()
+        torch.manual_seed(42)
+        derivative.simulate(n_paths=1)
+        with pytest.raises(ValueError):
+            m.implied_volatility()
+        result = m.implied_volatility(price=derivative.underlier.spot)
+        expect = m2.implied_volatility(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.spot,
+        )
+        assert_close(result, expect)
+        result = m.implied_volatility(
+            None,
+            derivative.max_log_moneyness(),
+            derivative.time_to_maturity(),
+            derivative.underlier.spot,
+        )
+        assert_close(result, expect)
+        result = m.implied_volatility(
+            derivative.log_moneyness(),
+            None,
+            derivative.time_to_maturity(),
+            derivative.underlier.spot,
+        )
+        assert_close(result, expect)
+        result = m.implied_volatility(
+            derivative.log_moneyness(),
+            derivative.max_log_moneyness(),
+            None,
+            derivative.underlier.spot,
+        )
+        assert_close(result, expect)
+        with pytest.raises(ValueError):
+            m2.implied_volatility(
+                None,
+                derivative.max_log_moneyness(),
+                derivative.time_to_maturity(),
+                derivative.underlier.spot,
+            )
+        with pytest.raises(ValueError):
+            m2.implied_volatility(
+                derivative.log_moneyness(),
+                None,
+                derivative.time_to_maturity(),
+                derivative.underlier.spot,
+            )
+        with pytest.raises(ValueError):
+            m2.implied_volatility(
+                derivative.log_moneyness(),
+                derivative.max_log_moneyness(),
+                None,
+                derivative.underlier.spot,
+            )
+
     def test_vega_and_gamma(self):
         m = BSLookbackOption()
         # vega = spot^2 * sigma * (T - t) * gamma
@@ -262,6 +824,27 @@ def test_vega_and_gamma(self):
         gamma = m.gamma(spot.log(), spot.log(), t, v)
         assert_close(vega, spot.square() * v * t * gamma, atol=1e-3, rtol=0)
 
+    @pytest.mark.parametrize("call", [True])
+    def test_vega_and_gamma_2(self, call: bool):
+        derivative = LookbackOption(BrownianStock(), call=call)
+        m = BSLookbackOption.from_derivative(derivative)
+        torch.manual_seed(42)
+        derivative.simulate(n_paths=1)
+        vega = m.vega()
+        gamma = m.gamma()
+        # ToDo: [..., :-1] should be removed
+        assert_close(
+            vega[..., :-1],
+            (
+                derivative.underlier.spot.square()
+                * derivative.underlier.volatility
+                * derivative.time_to_maturity()
+                * gamma
+            )[..., :-1],
+            atol=1e-3,
+            rtol=0,
+        )
+
     def test_put_notimplemented(self):
         with pytest.raises(ValueError):
             # not yet supported
diff --git a/tests/nn/test_functional.py b/tests/nn/test_functional.py
index bb3f5876..357ee647 100644
--- a/tests/nn/test_functional.py
+++ b/tests/nn/test_functional.py
@@ -4,6 +4,7 @@
 import torch
 from torch.testing import assert_close
 
+from pfhedge.nn.functional import bilerp
 from pfhedge.nn.functional import clamp
 from pfhedge.nn.functional import d1
 from pfhedge.nn.functional import d2
@@ -293,3 +294,31 @@ def test_d2_2():
             time_to_maturity=torch.as_tensor([1.0, 0, 0]),
             volatility=torch.as_tensor([0.2, 0.2, -1.0]),
         )
+
+
+def test_bilerp():
+    torch.manual_seed(42)
+
+    i1 = torch.randn(2, 3)
+    i2 = torch.randn(2, 3)
+    i3 = torch.randn(2, 3)
+    i4 = torch.randn(2, 3)
+
+    # edge cases
+    result = bilerp(i1, i2, i3, i4, 0.0, 0.0)
+    assert_close(result, i1)
+    result = bilerp(i1, i2, i3, i4, 1.0, 0.0)
+    assert_close(result, i2)
+    result = bilerp(i1, i2, i3, i4, 0.0, 1.0)
+    assert_close(result, i3)
+    result = bilerp(i1, i2, i3, i4, 1.0, 1.0)
+    assert_close(result, i4)
+
+    # w1 or w2 = 0 reduces to lerp
+    result = bilerp(i1, i2, i3, i4, 0.1, 0.0)
+    assert_close(result, torch.lerp(i1, i2, 0.1))
+    result = bilerp(i1, i2, i3, i4, 0.0, 0.1)
+    assert_close(result, torch.lerp(i1, i3, 0.1))
+
+    result = bilerp(i1, i2, i3, i4, 0.5, 0.5)
+    assert_close(result, (i1 + i2 + i3 + i4) / 4)
diff --git a/tests/test_version.py b/tests/test_version.py
new file mode 100644
index 00000000..8eda65cf
--- /dev/null
+++ b/tests/test_version.py
@@ -0,0 +1,10 @@
+import configparser
+import os.path
+
+import pfhedge
+
+
+def test_version():
+    parser = configparser.ConfigParser()
+    parser.read(os.path.join(os.path.dirname(__file__), "..", "pyproject.toml"))
+    assert pfhedge.__version__ == parser["tool.poetry"]["version"].replace('"', "")