diff --git a/vega_query/service/utils/party.py b/vega_query/service/utils/party.py index f1bb3d55..3db82514 100644 --- a/vega_query/service/utils/party.py +++ b/vega_query/service/utils/party.py @@ -60,6 +60,7 @@ def historic_balances( self, party_id: str, asset_id: str, + asset_decimals: int, market_id: Optional[str] = None, account_types: Optional[List[protos.vega.vega.AccountType.Value]] = None, date_range_start_timestamp: Optional[int] = None, @@ -102,4 +103,6 @@ def historic_balances( df = pd.DataFrame.from_dict(data, orient="index").sort_index().ffill() df["total"] = df.sum(axis=1) + df["total"] = df["total"]* 10**-asset_decimals return df + diff --git a/vega_query/visualisations/overlay.py b/vega_query/visualisations/overlay.py index 9b299454..15c07a41 100644 --- a/vega_query/visualisations/overlay.py +++ b/vega_query/visualisations/overlay.py @@ -30,7 +30,7 @@ def overlay_mark_price( ax: Axes, market_data_history: List[protos.vega.vega.MarketData], price_decimals: int, -): +)-> pd.DataFrame: x = [] y = [] for market_data in market_data_history: @@ -39,13 +39,19 @@ def overlay_mark_price( y.append(price if price != 0 else np.nan) ax.step(x, y, label="mark_price", where="post") + # Create a DataFrame from the extracted data + df_mark_price = pd.DataFrame({'timestamp': x, 'last_traded_price': y}) + + # Return the DataFrame + return df_mark_price + def overlay_last_traded_price( ax: Axes, market_data_history: List[protos.vega.vega.MarketData], price_decimals: int, **kwargs, -): +) : x = [] y = [] for market_data in market_data_history: diff --git a/vega_query/visualisations/plots/amm.py b/vega_query/visualisations/plots/amm.py index 4fdddd5f..4301864e 100644 --- a/vega_query/visualisations/plots/amm.py +++ b/vega_query/visualisations/plots/amm.py @@ -69,7 +69,7 @@ def create( party_colors = __party_colors(party_ids) - fig = plt.figure(figsize=(15, 8)) + fig = plt.figure(figsize=(16, 9)) gs = fig.add_gridspec(2, 2) ymin = ymax = 0 axes: List[Axes] = [] @@ -77,14 +77,15 @@ def create( axn0l = fig.add_subplot(gs[0, 0]) axn0r: Axes = axn0l.twinx() if market_data_history is not None: - overlay_mark_price(axn0l, market_data_history, market.decimal_places) + df_markPrice=overlay_mark_price(axn0l, market_data_history, market.decimal_places) + df_markPrice.to_csv("mark_price_data.csv", index=False) + overlay_last_traded_price(axn0l, market_data_history, market.decimal_places) overlay_trading_mode(axn0r, market_data_history) - # overlay_auction_starts(axn0r, market_data_history) - # overlay_auction_ends(axn0r, market_data_history) axn0l.set_ylabel("USDT") axn0l.set_title( - "Market Price: Ornstein–Uhlenbeck ($\\theta=0.01$, $\\sigma=0.5$)", + # "Market Price: Ornstein–Uhlenbeck ($\\theta=0.01$, $\\sigma=0.5$)", + "Market Price: Ornstein–Uhlenbeck", loc="left", ) @@ -106,7 +107,7 @@ def create( axes.append(ax11) ax21 = fig.add_subplot(gs[1, 1]) ax21.set_title("AMM: Aggregated Account Balances [USDT]", loc="left") - ax21.set_ylabel("$\\Delta$ USDT") + ax21.set_ylabel("USDT") ax21.sharex(axn0l) axes.append(ax21) @@ -141,7 +142,8 @@ def create( party_id=amm_party_id, size_decimals=market.position_decimal_places, ) - ax11.get_lines()[-1].set_label(amm_party_id[:7]) + # ax11.get_lines()[-1].set_label(amm_party_id[:7]) + ax11.get_lines()[-1].set_label("AMM") # Reconstruct the AMMs aggregated balance from balance changes df = service.utils.party.historic_balances( @@ -149,8 +151,11 @@ def create( asset_id=asset.id, date_range_start_timestamp=start_timestamp, date_range_end_timestamp=end_timestamp, + asset_decimals=asset.details.decimals, ) - ax21.step(df.index, df.total, where="post", label="total") + ax21.step(df.index, df.total, where="post", label="AMM portfolio") + df.to_csv("AMM_Portfolio_data.csv", index=False) + ax11.legend() ax21.legend() diff --git a/vega_sim/scenario/amm/registry.py b/vega_sim/scenario/amm/registry.py index 14fb3043..fb043d93 100644 --- a/vega_sim/scenario/amm/registry.py +++ b/vega_sim/scenario/amm/registry.py @@ -10,11 +10,12 @@ benchmark_configs=[ BenchmarkConfig( market_config=configs.mainnet.BTCUSDT.CONFIG, - initial_price=70000, - annualised_volatility=0.5, - notional_trade_volume=100, - process_theta=0.01, - process_drift=-10, + initial_price=60000, + annualised_volatility=0.28, + notional_trade_volume=10, + process_theta=0.0001, + process_drift=-5, + risky_trader_funds=1, ), ], amm_liquidity_fee=0.0001, diff --git a/vega_sim/scenario/amm/scenario.py b/vega_sim/scenario/amm/scenario.py index 317a38ff..d6bc945a 100644 --- a/vega_sim/scenario/amm/scenario.py +++ b/vega_sim/scenario/amm/scenario.py @@ -9,10 +9,11 @@ from vega_sim.scenario.common.agents import ( StateAgent, AutomatedMarketMaker, + MarketOrderTrader, + LimitOrderTrader, ExponentialShapedMarketMaker, ) - class AMMScenario(BenchmarkScenario): def __init__( self, @@ -58,6 +59,22 @@ def configure_agents( # Set commitment amount to 0 to avoid benchmark market maker from # providing liquidity, they should only supply volume through orders. agent.fee_amount = self.amm_liquidity_fee + agent.commitment_amount = 0.000001 + agent.supplied_amount = 0 + + # if isinstance(agent, MarketOrderTrader): + # # Set commitment amount to 0 to avoid benchmark market maker from + # # providing liquidity, they should only supply volume through orders. + # agent.fee_amount = self.amm_liquidity_fee + # agent.commitment_amount =10 + # agent.supplied_amount = 50 + + # if isinstance(agent, LimitOrderTrader): + # # Set commitment amount to 0 to avoid benchmark market maker from + # # providing liquidity, they should only supply volume through orders. + # agent.spread = 0 + # agent.mean = -7.5 + # agent.sigma = 0.3 # For each market, add an AMM agent. for benchmark_config in self.benchmark_configs: @@ -67,8 +84,8 @@ def configure_agents( wallet_name="AutomatedMarketMaker", key_name=f"AutomatedMarketMaker_{benchmark_config.market_config.instrument.code}_{str(i_agent).zfill(3)}", market_name=benchmark_config.market_config.instrument.name, - initial_asset_mint=1e6, - commitment_amount=7000, + initial_asset_mint=1e5, + commitment_amount=6000, slippage_tolerance=0.05, proposed_fee=self.amm_liquidity_fee, price_process=iter(benchmark_config.price_process), diff --git a/vega_sim/scenario/ammV2/registryV2.py b/vega_sim/scenario/ammV2/registryV2.py new file mode 100644 index 00000000..57646acd --- /dev/null +++ b/vega_sim/scenario/ammV2/registryV2.py @@ -0,0 +1,28 @@ +import vega_sim.configs as configs +from vega_sim.scenario.benchmark.configs import BenchmarkConfig +from vega_sim.scenario.ammV2.scenarioV2 import AMMScenarioV2 + + +REGISTRYV2 = { + "static": AMMScenarioV2( + block_length_seconds=1, + step_length_seconds=30, + benchmark_configs=[ + BenchmarkConfig( + market_config=configs.mainnet.BTCUSDT.CONFIG, + initial_price=60000, + annualised_volatility=0.28, + notional_trade_volume=800, + process_theta=0.0001, + process_drift=-5, + risky_trader_funds=1, + ), + ], + amm_liquidity_fee=0.0001, + amm_update_frequency=0, + initial_network_parameters={ + "validators.epoch.length": "1h", + "market.fee.factors.makerFee": "0", + }, + ), +} diff --git a/vega_sim/scenario/ammV2/runV2.py b/vega_sim/scenario/ammV2/runV2.py new file mode 100644 index 00000000..4928ff8e --- /dev/null +++ b/vega_sim/scenario/ammV2/runV2.py @@ -0,0 +1,107 @@ +import os +import logging +import pathlib +import datetime +import argparse + + +from vega_sim.null_service import VegaServiceNull, Ports +from vega_sim.scenario.constants import Network +from vega_sim.scenario.ammV2.scenarioV2 import AMMScenarioV2 +from vega_sim.scenario.ammV2.registryV2 import REGISTRYV2 + +from vega_query.service.service import Service +from vega_query.service.networks.constants import Network + +import vega_query.visualisations as vis +from matplotlib.figure import Figure + + +def _run( + scenario: AMMScenarioV2, + pause: bool = False, + console: bool = False, + output: bool = False, + wallet: bool = False, + output_dir: str = "plots", + core_metrics_port: int = 2723, + data_node_metrics_port: int = 3651, +): + + with VegaServiceNull( + warn_on_raw_data_access=False, + seconds_per_block=scenario.block_length_seconds, + transactions_per_block=scenario.transactions_per_block, + retain_log_files=True, + use_full_vega_wallet=wallet, + run_with_console=console, + port_config={ + Ports.METRICS: core_metrics_port, + Ports.DATA_NODE_METRICS: data_node_metrics_port, + }, + ) as vega: + scenario.run_iteration( + vega=vega, + log_every_n_steps=100, + output_data=False, + run_with_snitch=False, + ) + + if output: + if not os.path.exists(output_dir): + os.mkdir(output_dir) + output_dir = output_dir + f"/{datetime.datetime.now()}" + if not os.path.exists(output_dir): + os.mkdir(output_dir) + + # Use vegapy package to produce plots + service = Service( + network=Network.NETWORK_LOCAL, + port_data_node=vega.data_node_grpc_port, + ) + + for benchmark_config in scenario.benchmark_configs: + fig: Figure = vis.plots.amm.create( + service, + market_code=benchmark_config.market_config.instrument.code, + ) + fig.savefig( + f"{output_dir}/{benchmark_config.market_config.instrument.code.replace('/','-')}_amm.png" + ) + + if pause: + input("Waiting after run finished.") + + +if __name__ == "__main__": + + parser = argparse.ArgumentParser() + parser.add_argument("-m", "--market", required=True, type=str) + parser.add_argument("-s", "--steps", default=600, type=int) + parser.add_argument("-p", "--pause", action="store_true") + parser.add_argument("-d", "--debug", action="store_true") + parser.add_argument("-o", "--output", action="store_true") + parser.add_argument("-c", "--console", action="store_true") + parser.add_argument("-w", "--wallet", action="store_true") + parser.add_argument("--core-metrics-port", default=2723, type=int) + parser.add_argument("--data-node-metrics-port", default=3651, type=int) + args = parser.parse_args() + + logging.basicConfig( + level=logging.DEBUG if args.debug else logging.INFO, + format="%(asctime)s - %(name)s - %(levelname)s - %(message)s", + ) + + if args.market not in REGISTRYV2: + raise ValueError(f"Market {args.market} not found") + scenario = REGISTRYV2[args.market].num_steps = args.steps + + _run( + scenario=REGISTRYV2[args.market], + wallet=args.wallet, + console=args.console, + pause=args.pause, + output=args.output, + core_metrics_port=args.core_metrics_port, + data_node_metrics_port=args.data_node_metrics_port, + ) diff --git a/vega_sim/scenario/ammV2/scenarioV2.py b/vega_sim/scenario/ammV2/scenarioV2.py new file mode 100644 index 00000000..97b18af6 --- /dev/null +++ b/vega_sim/scenario/ammV2/scenarioV2.py @@ -0,0 +1,83 @@ +import itertools +import numpy as np +import re +from typing import Optional, List, Dict, Any + +from vega_sim.null_service import VegaServiceNull +from vega_sim.scenario.benchmark.configs import BenchmarkConfig +from vega_sim.scenario.benchmark.scenario_Informed import BenchmarkScenario_Informed +from vega_sim.scenario.common.agents import ( + StateAgent, + AutomatedMarketMaker, + MarketOrderTrader, + LimitOrderTrader, + InformedTrader, +) + +class AMMScenarioV2(BenchmarkScenario_Informed): + def __init__( + self, + benchmark_configs: List[BenchmarkConfig], + num_steps: int = 60 * 24 * 30 * 3, + transactions_per_block: int = 4096, + block_length_seconds: float = 1, + step_length_seconds: Optional[float] = None, + amm_liquidity_fee: float = 0.0001, + amm_update_frequency: float = 0.1, + output: bool = True, + initial_network_parameters: Optional[Dict[str, Any]] = None, + ): + super().__init__( + benchmark_configs=benchmark_configs, + num_steps=num_steps, + transactions_per_block=transactions_per_block, + block_length_seconds=block_length_seconds, + step_length_seconds=step_length_seconds, + output=output, + initial_network_parameters=initial_network_parameters, + ) + + self.amm_liquidity_fee = amm_liquidity_fee + self.amm_update_frequency = amm_update_frequency + + def configure_agents( + self, + vega: VegaServiceNull, + tag: str, + random_state: Optional[np.random.RandomState], + **kwargs, + ) -> Dict[str, StateAgent]: + self.random_state = ( + random_state if random_state is not None else np.random.RandomState() + ) + agents = super().configure_agents(vega, tag, random_state, **kwargs) + extra_agents = [] + + agents = super().configure_agents(vega, tag, random_state, **kwargs) + + # For each market, add an AMM agent. + for benchmark_config in self.benchmark_configs: + i_agent = 0 + extra_agents.append( + AutomatedMarketMaker( + wallet_name="AutomatedMarketMaker", + key_name=f"AutomatedMarketMaker_{benchmark_config.market_config.instrument.code}_{str(i_agent).zfill(3)}", + market_name=benchmark_config.market_config.instrument.name, + initial_asset_mint=1e5, + commitment_amount=6000, + slippage_tolerance=0.05, + proposed_fee=self.amm_liquidity_fee, + price_process=iter(benchmark_config.price_process), + lower_bound_scaling=1 - 0.02, + upper_bound_scaling=1 + 0.02, + leverage_at_lower_bound=20, + leverage_at_upper_bound=20, + update_bias=self.amm_update_frequency, + tag=f"{benchmark_config.market_config.instrument.code}_{str(i_agent).zfill(3)}", + random_state=self.random_state, + ) + ) + + extra_agents = {agent.name(): agent for agent in extra_agents} + agents.update(extra_agents) + return agents diff --git a/vega_sim/scenario/benchmark/scenario.py b/vega_sim/scenario/benchmark/scenario.py index ffb33c6c..bc97d2c0 100644 --- a/vega_sim/scenario/benchmark/scenario.py +++ b/vega_sim/scenario/benchmark/scenario.py @@ -154,8 +154,7 @@ def configure_agents( buy_intensity=100, sell_intensity=100, base_order_size=benchmark_config.notional_trade_volume - / benchmark_config.initial_price - / 100, + / (benchmark_config.initial_price*100), step_bias=1, initial_asset_mint=1e8, tag=f"{benchmark_config.market_config.instrument.code}_{str(i_agent).zfill(3)}", @@ -171,11 +170,9 @@ def configure_agents( market_name=market_name, time_in_force_opts={"TIME_IN_FORCE_GTT": 1}, buy_volume=benchmark_config.notional_trade_volume - / benchmark_config.initial_price - / 100, + / (benchmark_config.initial_price*100), sell_volume=benchmark_config.notional_trade_volume - / benchmark_config.initial_price - / 100, + / (benchmark_config.initial_price*100), buy_intensity=100, sell_intensity=100, submit_bias=1, diff --git a/vega_sim/scenario/benchmark/scenario_Informed.py b/vega_sim/scenario/benchmark/scenario_Informed.py new file mode 100644 index 00000000..b41ba7b9 --- /dev/null +++ b/vega_sim/scenario/benchmark/scenario_Informed.py @@ -0,0 +1,236 @@ +import numpy as np +from typing import Optional, Dict, Any, List + +from vega_sim.api.market import MarketConfig +from vega_sim.scenario.scenario import Scenario +from vega_sim.scenario.benchmark.configs import BenchmarkConfig +from vega_sim.scenario.common.utils.price_process import ou_price_process +from vega_sim.scenario.constants import Network +from vega_sim.null_service import VegaServiceNull +from vega_sim.environment.environment import ( + MarketEnvironmentWithState, + NetworkEnvironment, +) + +from vega_sim.configs.agents import ConfigurableMarketManager +from vega_sim.scenario.common.agents import ( + StateAgent, + NetworkParameterManager, + UncrossAuctionAgent, + InformedTrader, + MarketOrderTrader, + LimitOrderTrader, +) +from vega_sim.scenario.fuzzed_markets.agents import ( + RiskyMarketOrderTrader, +) + + +class BenchmarkScenario_Informed(Scenario): + def __init__( + self, + benchmark_configs: List[BenchmarkConfig], + num_steps: int = 60 * 24 * 30 * 3, + transactions_per_block: int = 4096, + block_length_seconds: float = 1, + step_length_seconds: Optional[float] = None, + initial_network_parameters: Dict[str, Any] = None, + output: bool = True, + ): + super().__init__() + + self.num_steps = num_steps + self.step_length_seconds = ( + step_length_seconds + if step_length_seconds is not None + else block_length_seconds + ) + self.block_length_seconds = block_length_seconds + self.transactions_per_block = transactions_per_block + + self.output = output + self.benchmark_configs = benchmark_configs + + self.initial_network_parameters = ( + initial_network_parameters if initial_network_parameters is not None else {} + ) + + def configure_agents( + self, + vega: VegaServiceNull, + tag: str, + random_state: Optional[np.random.RandomState], + **kwargs, + ) -> Dict[str, StateAgent]: + self.random_state = ( + random_state if random_state is not None else np.random.RandomState() + ) + + self.agents = [] + self.agents.append( + NetworkParameterManager( + wallet_name="NetworkParameterManager", + key_name="NetworkParameterManager", + network_parameters=self.initial_network_parameters, + ) + ) + for _, benchmark_config in enumerate(self.benchmark_configs): + + market_name = benchmark_config.market_config.instrument.name + market_decimal_places = int( + benchmark_config.market_config.price_decimal_places + if benchmark_config.market_config.is_spot() + else benchmark_config.market_config.decimal_places + ) + # Create fuzzed price process + benchmark_config.price_process = ou_price_process( + self.num_steps + 1, + x0=benchmark_config.initial_price, + mu=benchmark_config.initial_price, + theta=benchmark_config.process_theta, + drift=benchmark_config.process_drift, + sigma=benchmark_config.annualised_volatility + * np.sqrt(self.step_length_seconds / (365.25 * 24 * 60 * 60)) + * benchmark_config.initial_price, + ) + self.agents.append( + ConfigurableMarketManager( + wallet_name="ConfigurableMarketManager", + key_name=f"ConfigurableMarketManager_{benchmark_config.market_config.instrument.code}", + market_config=benchmark_config.market_config, + oracle_prices=iter(benchmark_config.price_process), + oracle_submission=0.5, + oracle_difference=0.001, + random_state=self.random_state, + tag=f"{benchmark_config.market_config.instrument.code}", + ), + ) + + self.agents.append( + InformedTrader( + wallet_name="InformedTrader", + key_name=f"InformedTrader_{benchmark_config.market_config.instrument.code}", + # price_process_generator=iter(benchmark_config.price_process), + price_process=benchmark_config.price_process, + initial_asset_mint=1e9, + market_name=market_name, + proportion_taken=0.8, + accuracy=1, + lookahead=5, + max_abs_position=100, + tag=f"{benchmark_config.market_config.instrument.code}", + ) + ) + + self.agents.extend( + [ + UncrossAuctionAgent( + wallet_name="UncrossAuctionAgent", + key_name=f"{benchmark_config.market_config.instrument.code}_{side}", + side=side, + initial_asset_mint=1e8, + price_process=iter(benchmark_config.price_process), + market_name=market_name, + uncrossing_size=1e6 / benchmark_config.initial_price, + tag=( + f"{benchmark_config.market_config.instrument.code}_{str(i_agent).zfill(3)}" + ), + leave_opening_auction_prob=1.0, + ) + for i_agent, side in enumerate(["SIDE_BUY", "SIDE_SELL"]) + ] + ) + self.agents.extend( + [ + MarketOrderTrader( + wallet_name="MarketOrderTrader", + key_name=f"MarketOrderTrader_{str(i_agent).zfill(3)}", + market_name=market_name, + buy_intensity=100, + sell_intensity=100, + base_order_size=benchmark_config.notional_trade_volume + / (benchmark_config.initial_price*100), + step_bias=1, + initial_asset_mint=1e8, + tag=f"{benchmark_config.market_config.instrument.code}_{str(i_agent).zfill(3)}", + ) + for i_agent in range(1) + ] + ) + self.agents.extend( + [ + LimitOrderTrader( + wallet_name=f"LimitOrderTrader", + key_name=f"LimitOrderTrader_{str(i_agent).zfill(3)}", + market_name=market_name, + time_in_force_opts={"TIME_IN_FORCE_GTT": 1}, + buy_volume=benchmark_config.notional_trade_volume + / (benchmark_config.initial_price*100), + sell_volume=benchmark_config.notional_trade_volume + / (benchmark_config.initial_price*100), + buy_intensity=100, + sell_intensity=100, + submit_bias=1, + cancel_bias=0, + duration=120, + price_process=benchmark_config.price_process, + spread=0, + mean=-3, + sigma=0.5, + initial_asset_mint=1e8, + tag=f"{benchmark_config.market_config.instrument.code}_{str(i_agent).zfill(3)}", + ) + for i_agent in range(1) + ] + ) + + if ( + benchmark_config.market_config.is_future() + or benchmark_config.market_config.is_perp() + ): + self.agents.extend( + [ + RiskyMarketOrderTrader( + wallet_name="RiskyMarketOrderTrader", + key_name=f"RiskyMarketOrderTrader_{benchmark_config.market_config.instrument.code}_{str(i_agent).zfill(3)}", + market_name=market_name, + side=side, + initial_asset_mint=benchmark_config.risky_trader_funds, + leverage_factor=0.5, + step_bias=0.1, + tag=f"{benchmark_config.market_config.instrument.code}_{side}_{str(i_agent).zfill(3)}", + ) + for side in ["SIDE_BUY", "SIDE_SELL"] + for i_agent in range(5) + ] + ) + + return {agent.name(): agent for agent in self.agents} + + def configure_environment( + self, + vega: VegaServiceNull, + **kwargs, + ) -> MarketEnvironmentWithState: + if kwargs.get("network", Network.NULLCHAIN) == Network.NULLCHAIN: + return MarketEnvironmentWithState( + agents=list(self.agents.values()), + n_steps=self.num_steps, + random_agent_ordering=False, + transactions_per_block=self.transactions_per_block, + vega_service=vega, + step_length_seconds=self.step_length_seconds, + block_length_seconds=vega.seconds_per_block, + ) + else: + return NetworkEnvironment( + agents=list(self.agents.values()), + n_steps=self.num_steps, + vega_service=vega, + step_length_seconds=self.step_length_seconds, + raise_datanode_errors=kwargs.get("raise_datanode_errors", False), + raise_step_errors=kwargs.get("raise_step_errors", False), + random_state=self.random_state, + create_keys=True, + mint_keys=True, + ) diff --git a/vega_sim/scenario/common/utils/price_process.py b/vega_sim/scenario/common/utils/price_process.py index 39e72b29..aec39437 100644 --- a/vega_sim/scenario/common/utils/price_process.py +++ b/vega_sim/scenario/common/utils/price_process.py @@ -62,7 +62,7 @@ def random_walk( S = np.round(S, decimal_precision) # If random state is passed then error if it generates a negative price - # Otherwise retry with a new seed + # Otherwise retry with a new if trim_to_min is not None: S[S < trim_to_min] = trim_to_min @@ -315,7 +315,7 @@ def get_live_price( return _live_prices[feed_key] -def ou_price_process(n, theta=0.15, mu=0.0, sigma=0.2, x0=1.0, drift=0.0): +def ou_price_process(n, theta=0.15, mu=0.0, sigma=0.2, x0=1.0, drift=0.0, seed=42): """ Generates a mean-reverting price series using the Ornstein–Uhlenbeck process with an optional drift term. @@ -331,6 +331,8 @@ def ou_price_process(n, theta=0.15, mu=0.0, sigma=0.2, x0=1.0, drift=0.0): Returns: np.ndarray: Mean-reverting price series of length n. """ + if seed is not None: + np.random.seed(seed) # Set the random seed for reproducibility dt = 1.0 # Assuming a time step of 1 x = np.zeros(n) x[0] = x0