regression parameter estimation for correlated variables Usually we don't want to include correlated variables in regression model as problems with estimation and variable significance arise.
I have always thought that a major problem is that the estimates of parameters of correlated variables will "split the total effect". And so the estimates will be (very) different from those obtained if only one of the correlated variables is included in the model.
However, when I simulated a situation where two almost perfectly correlated variables enter the regression model, it turned out that most of the times the true parameter is estimated very well. So, my question is - how the algorithm knows which is the true parameter to estimate, if both variables are very (!) similar, and why does the effect is not "split" between the variables?
library(MASS)
library(tidyverse)

rep_res<-replicate(1000,{
r <- 0.99 # almost perfectly correlated variables
samples <- 200
var_df <-  mvrnorm(n=samples, mu=c(0, 0),
                   Sigma=matrix(c(1, r, r, 1), nrow=2)) %>% as.data.frame()
# Y depends only on V1
# true regression equation y = 2 + 3*V1
Y <- 2 + 3*var_df$V1 + rnorm(samples,0,1)
res <- lm(Y ~ V2 + V1, data = var_df) # include V1 and V2 in model
res %>% summary() %>% coef() %>% .[2:3,1]})

res_df <- t(rep_res) %>% as.data.frame() %>% melt()

Below is the coefficients from 1000 replications. 

 A: 
Usually we don't want to include correlated variables in regression model as problems with estimation and variable significance arise.

This is not true. The situation in which all regressors are uncorrelated is a special case, multiple regression is a powerful tool also because does not require so strict condition. Is right that correlations among regressors inflate variance of estimators. Sample size matters but rarely mid level correlations take away any value at the estimators.

I have always thought that a major problem is that the estimates of parameters of correlated variables will "split the total effect". And so the estimates will be (very) different from those obtained if only one of the correlated variables is included in the model.

You have to keep in mind the what you call "total effect" is the usual correlation. Is true that, in certain sense, multiple regression permit us to split the correlation/effect but this is precisely what we want.
The case of perfect collinearity, then perfect correlation (linear), is a limit situation. If it do not appear, the well known OLS formula work:
$\beta_{OLS} = (X'X)^{-1} X'y$
If the true model contain precisely the variables collected in $X$, or more, the OLS estimator converge to the true parameters. If $X$ contain less variables you encounter the omitted variables problem.
A: You variables are not perfectly correlated (r=1) as you stated.
Suppose that $V_1$ and $V_2$ are two correlated , zero mean, unit variance random variables. 
$$E(V_1)=E(V_2)=0$$
$$E(V^2_1)=E(V^2_2)=1$$
$$E(V_1 V_2)=r$$
The true model of the data is
$$X=a_1 V_1+ n$$
where n is an additive noise component.
you are however trying to model the data as
$$\hat X=\hat a_1 V_1 +\hat a_2 V_2 $$.
If you solve the mean square estimation problem by minimizing 
$$\frac{d}{d \hat a_1} E(a_1 v_1 - \hat a_1 V_1 -\hat a_2 V_2)^2 =0$$
$$\frac{d}{d \hat a_2} E(a_1 v_1 - \hat a_1 V_1 -\hat a_2 V_2)^2 =0$$
you obtain the following system of equations:
$$\hat a_1  +\hat a_2 r = a_1$$
$$\hat a_1 r +\hat a_2  = a_1 r$$
which leads to 
$\hat a_1 = a_1 \frac{1-r^2}{1-r^2}$.
So as long as you don't have perfect correlation you can in principle recover the correct coefficients. 
One practical problem is that you are actually solving a least square problem where the expected values are replaced by summation. In this case projection of $V_1$ and $V_2$ along the noise component $n$
$$\sum^N_{i=1} {V_{1i} n_i}$$
$$\sum^N_{i=1} {V_{2i} n_i}$$
would be non zero and will be amplified when $r$ is close to 1 (owing to division by $1-r^2$).
