# set parameters
	X1_shares <- 10
	X2_shares <- -2
	P1_start <- 2
	P2_start <- 2
	T <- 100  # number of time periods in each trial
	numb_trials <- 100000
	set.seed(12345)
	
# create covariance matrix
	Sigma <- matrix(NA,2,2)
	Sigma[1,1] <- 0.03863539
	Sigma[2,2] <- 0.05585887
	Sigma[1,2] <- 0.02213511
	Sigma[2,1] <- Sigma[1,2]
	cat("\nSigma:\n")
	print(Sigma)

# function that does one trial
do_trial <- function(T,X1_shares,X2_shares,P1_start,P2_start,Sigma) {
	# get choleski decomp of covariance matrix
		U <- chol(Sigma)
	# set up container to hold values of P1 and P2 over all time periods
	# one row for each time period, one column for each variable
		Price_matrix <- matrix(NA,T,2)
	# initial values for P1 and P2
		P1 <- P1_start
		P2 <- P2_start
	# loop through time periods to build up P1 and P2 series
		for (t in 0:T) {
			# generate correlated returns
			returns <- U%*%rnorm(2)
			# generate new values for P1 and P2
			P1 <- P1*exp(returns[1])
			P2 <- P2*exp(returns[2])
			# save new values in container
			Price_matrix[t,1] <- P1
			Price_matrix[t,2] <- P2
			}
	# loop through values, check for ruin
		ruin <- FALSE
		for (t in 1:T) {
			if (ruin == FALSE) {
				if (X1_shares*Price_matrix[t,1]+X2_shares*Price_matrix[t,2] < 0) {
					ruin <- TRUE
				}
			}
		}
	return(ruin)
	}

# function that does one trial version 2
do_trial_v2 <- function(T,X1_shares,X2_shares,P1_start,P2_start,Sigma) {
	# get choleski decomp of covariance matrix
		U <- chol(Sigma)
	# set up container to hold values of P1 and P2 over all time periods
	# one row for each time period, one column for each variable
		Price_matrix <- matrix(NA,T,2)
	# initial values for P1 and P2
		P1 <- P1_start
		P2 <- P2_start
	# create P1 and P2 series using single call to rnorm and no loop
		# draw all 2*T deviates simultaneously
		ISND <- rnorm(2*T)
		# dimension correctly for matric multiplecation
		dim(ISND) <- c(2,T)
		# correlate them (U is 2x2, ISDN is 2xT, CSND is 2xT; one row for each price series)
		CSND <- U%*%ISND 
		# cummulative products
		P1_cummulative_products <- cumprod(exp(CSND[1,]))	
		P2_cummulative_products <- cumprod(exp(CSND[2,]))	
		Price_matrix[,1] <- P1_start*P1_cummulative_products
		Price_matrix[,2] <- P2_start*P2_cummulative_products
	# loop through values, check for ruin
		ruin <- FALSE
		for (t in 1:T) {
			if (ruin == FALSE) {
				if (X1_shares*Price_matrix[t,1]+X2_shares*Price_matrix[t,2] < 0) {
					ruin <- TRUE
				}
			}
		}
	return(ruin)
	}

# function that does one trial version 3
do_trial_v3 <- function(T,X1_shares,X2_shares,P1_start,P2_start,Sigma) {
	# get choleski decomp of covariance matrix
		U <- chol(Sigma)
	# set up container to hold values of P1 and P2 over all time periods
	# one row for each time period, one column for each variable
		Price_matrix <- matrix(NA,T,2)
	# initial values for P1 and P2
		P1 <- P1_start
		P2 <- P2_start
	# create P1 and P2 series using single call to rnorm and no loop
		# draw all 2*T deviates simultaneously
		ISND <- rnorm(2*T)
		# dimension correctly for matric multiplecation
		dim(ISND) <- c(2,T)
		# correlate them (U is 2x2, ISDN is 2xT, CSND is 2xT; one row for each price series)
		CSND <- U%*%ISND 
		# cummulative products
		P1_cummulative_products <- cumprod(exp(CSND[1,]))	
		P2_cummulative_products <- cumprod(exp(CSND[2,]))	
		Price_matrix[,1] <- P1_start*P1_cummulative_products
		Price_matrix[,2] <- P2_start*P2_cummulative_products
	# loop through values, check for ruin
		ruin <- FALSE
		X <- matrix(NA,2,1)
		X[1,1] <- X1_shares
		X[2,1] <- X2_shares
		portfolio_value_vec <- Price_matrix %*% X
		if (min(portfolio_value_vec) < 0) ruin <- TRUE
	return(ruin)
	}



# do trials
	Rprof("profile_out.txt")

	ruin_sum <- 0
	for (i in 1:numb_trials) {
		#ruin_sum <- ruin_sum + do_trial(T,X1_shares,X2_shares,P1_start,P2_start,Sigma)
		#ruin_sum <- ruin_sum + do_trial_v2(T,X1_shares,X2_shares,P1_start,P2_start,Sigma)
		ruin_sum <- ruin_sum + do_trial_v3(T,X1_shares,X2_shares,P1_start,P2_start,Sigma)
	}
	cat("\nProbability of ruin:\n")
	print(ruin_sum/numb_trials)

	Rprof(NULL)
	summaryRprof("profile_out.txt")$by.self