From comment:

How the result of R is correct, because if you look at my run-though
of a program call, the value of Xs isn't [1,3], which is what it
eventually outputs; it is instead [3] which unifies to R (clearly I am
missing something along the way, but I am unsure what that is).

This is correct:

% Initial call

val_and_remainder(2, [1,2,3], R).



val_and_remainder(2, [1|[2,3]], [1|R]) :- val_and_remainder(2, [2,3], R).



% Hits base case

val_and_remainder(2, [2|[3]], [3]).

however Prolog is not like other programming languages where you enter with input and exit with output at a return statement. In Prolog you move forward through the predicate statements unifying and continuing with predicates that are true, and upon backtracking also unifying the unbound variables. (That is not technically correct but it is easier to understand for some if you think of it that way.)

You did not take into consideration the the unbound variables that are now bound upon backtracking.

When you hit the base case Xs was bound to [3],

but when you backtrack you have look at

val_and_remainder(2, [1|[2,3]], [1|R])

and in particular [1|R] for the third parameter.

Since Xs was unified with R in the call to the base case, i.e.

val_and_remainder(X,[X|Xs],Xs).

R now has [3].

Now the third parameter position in

val_and_remainder(2, [1|[2,3]], [1|R])

is [1|R] which is [1|[3]] which as syntactic sugar is [1,3] and not just [3].

Now when the query

val_and_remainder(2, [1,2,3], R).

was run, the third parameter of the query R was unified with the third parameter of the predicate

val_and_remainder(X,[Y|Ys],[Y|R])

so R was unified with [Y|R] which unpon backtracking is [1,3]
and thus the value bound to the query variable R is [1,3]

Stepwise reproduction of Prolog's mechanism often leads to more confusion than it helps. You probably have notions like "returning" meaning something very specific—more appropriate to imperative languages.

Here are different approaches you can always use:

Ask the most general query

... and let Prolog explain you what the relation is about.

| ?- val_and_remainder(X, Xs, Ys).

     Xs = [X|Ys]

  ;  Xs = [_A,X|_B],

     Ys = [_A|_B]

  ;  Xs = [_A,_B,X|_C],

     Ys = [_A,_B|_C]

  ;  Xs = [_A,_B,_C,X|_D],

     Ys = [_A,_B,_C|_D]

  ;  Xs = [_A,_B,_C,_D,X|_E],

     Ys = [_A,_B,_C,_D|_E]

  ; ...

So Xs and Ys share a common list prefix, Xs has thereafter an X, followed by a common rest. This query would continue producing further answers. Sometimes, you want to see all answers, then you have to be more specific. But don't be too specific:

| ?- Xs = [_,_,_,_], val_and_remainder(X, Xs, Ys).

     Xs = [X,_A,_B,_C],

     Ys = [_A,_B,_C]

  ;  Xs = [_A,X,_B,_C],

     Ys = [_A,_B,_C]

  ;  Xs = [_A,_B,X,_C],

     Ys = [_A,_B,_C]

  ;  Xs = [_A,_B,_C,X],

     Ys = [_A,_B,_C]

  ;  false.

So here we got all possible answers for a four-element list. All of them.

Stick to ground goals when going through specific inferences

So instead of val_and_remainder(2, [1,2,3], R). (which obviously got your head spinning) rather consider val_and_remainder(2, [1,2,3], [1,3]). and then
val_and_remainder(2, [2,3],[3]). From this side it should be obvious.

Read Prolog rules right-to-left

See Prolog rules as production rules. Thus, whenever everything holds on the right-hand side of a rule, you can conclude what is on the left. Thus, the :- is an early 1970s' representation of a ←

Later on, you may want to ponder more complex questions, too. Like

Functional dependencies

Does the first and second argument uniquely determine the last one? Does X, Xs → Ys hold?

Here is a sample query that asks for Ys and Ys2 being different for the same X and Xs.

| ?- val_and_remainder(X, Xs, Ys), val_and_remainder(X, Xs, Ys2), dif(Ys,Ys2).

     Xs = [X,_A,X|_B],

     Ys = [_A,X|_B],

     Ys2 = [X,_A|_B],

     dif([_A,X|_B],[X,_A|_B])

  ;  ...

So apparently, there are different values for Ys for a given X and Xs. Here is a concrete instance:

| ?- val_and_remainder(x, [x,a,x], Ys).

     Ys = [a,x]

  ;  Ys = [x,a]

  ;  false.

There is no classical returning here. It does not return once but twice. It's more of a yield.

Yet, there is in fact a functional dependency between the arguments! Can you find it? And can you Prolog-wise prove it (as much as Prolog can do a proof, indeed).

I don't understand the name of your predicate. It is a distraction anyway. The non-uniform naming of the variables is a distraction as well. Let's use some neutral, short one-syllable names to focus on the code itself in its clearest form:

foo( H, [H | T], T).                          % 1st clause



foo( X, [H | T], [H | R]) :- foo( X, T, R).   % 2nd clause

So it's the built-in select/3. Yay!..

Now you ask about the query foo( 2, [1,2,3], R) and how does R gets its value set correctly. The main thing missing from your rundown is the renaming of variables when a matching clause is selected. The resolution of the query goes like this:

|-  foo( 2, [1,2,3], R) ? { } 

  %% SELECT -- 1st clause, with rename

  |-  ?  { foo( H1, [H1|T1], T1) = foo( 2, [1,2,3], R) }

         **FAIL** (2 = 1)

         **BACKTRACK to the last SELECT**

  %% SELECT -- 2nd clause, with rename

  |-  foo( X1, T1, R1) ?

         { foo( X1, [H1|T1], [H1|R1]) = foo( 2, [1,2,3], R) }

         **OK**

    %% REWRITE

    |-  foo( X1, T1, R1) ?

          { X1=2, [H1|T1]=[1,2,3], [H1|R1]=R }

    %% REWRITE

    |-  foo( 2, [2,3], R1) ?  { R=[1|R1] } 

      %% SELECT -- 1st clause, with rename

      |-  ? { foo( H2, [H2|T2], T2) = foo( 2, [2,3], R1), R=[1|R1] }

             ** OK ** 

        %% REWRITE

        |-  ? { H2=2, T2=[3], T2=R1, R=[1|R1] }

        %% REWRITE

        |-  ? { R=[1,3] }

        %% DONE

The goals between |- and ? are the resolvent, the equations inside { } are the substitution. The knowledge base (KB) is implicitly to the left of |- in its entirety.

On each step, the left-most goal in the resolvent is chosen, a clause with the matching head is chosen among the ones in the KB ^{(while renaming all of the clause's variables in the consistent manner, such that no variable in the resolvent is used by the renamed clause, so there's no accidental variable capture)}, and the chosen goal is replaced in the resolvent with that clause's body, while the successful unification is added into the substitution. When the resolvent is empty, the query has been proven and what we see is the one successful and-branch in the whole and-or tree.

This is how a machine could be doing it. The "rewrite" steps are introduced here for ease of human comprehension.

So we can see here that the first successful clause selection results in the equation

     R = [1 | R1     ]

, and the second, --

              R1 = [3]

, which together entail

     R = [1,        3]

This gradual top-down instantiation / fleshing-out of lists is a very characteristic Prolog's way of doing things.

In response to the bounty challenge, regarding functional dependency in the relation foo/3 (i.e. select/3): in foo(A,B,C), any two ground values for B and C uniquely determine the value of A (or its absence):

2 ?- foo( A, [0,1,2,1,3], [0,2,1,3]).

A = 1 ;

false.



3 ?- foo( A, [0,1,2,1,3], [0,1,2,3]).

A = 1 ;

false.



4 ?- foo( A, [0,1,2,1,3], [0,1,2,4]).

false.



f ?- foo( A, [0,1,1], [0,1]).

A = 1 ;

A = 1 ;

false.

Attempt to disprove it by a counterargument:

10 ?- dif(A1,A2), foo(A1,B,C), foo(A2,B,C).

Action (h for help) ? abort

% Execution Aborted

Prolog fails to find a counterargument.

Tying to see more closely what's going on, with iterative deepening:

28 ?- length(BB,NN), foo(AA,BB,CC), XX=[AA,BB,CC], numbervars(XX), 

      writeln(XX), (NN>3, !, fail).

[A,[A],]

[A,[A,B],[B]]

[A,[B,A],[B]]

[A,[A,B,C],[B,C]]

[A,[B,A,C],[B,C]]

[A,[B,C,A],[B,C]]

[A,[A,B,C,D],[B,C,D]]

false.



29 ?- length(BB,NN), foo(AA,BB,CC), foo(AA2,BB,CC), 

      XX=[AA,AA2,BB,CC],  numbervars(XX), writeln(XX), (NN>3, !, fail).

[A,A,[A],]

[A,A,[A,B],[B]]

[A,A,[A,A],[A]]

[A,A,[A,A],[A]]

[A,A,[B,A],[B]]

[A,A,[A,B,C],[B,C]]

[A,A,[A,A,B],[A,B]]

[A,A,[A,A,A],[A,A]]

[A,A,[A,A,B],[A,B]]

[A,A,[B,A,C],[B,C]]

[A,A,[B,A,A],[B,A]]

[A,A,[A,A,A],[A,A]]

[A,A,[B,A,A],[B,A]]

[A,A,[B,C,A],[B,C]]

[A,A,[A,B,C,D],[B,C,D]]

false.

AA and AA2 are always instantiated to the same variable.

There's nothing special about the number 3, so it is safe to conjecture by generalization that it will always be so, for any length tried.

Another attempt at Prolog-wise proof:

ground_list(LEN,L):-

  findall(N, between(1,LEN,N), NS), 

  member(N,NS),

  length(L,N), 

  maplist( A^member(A,NS), L).



bcs(N, BCS):- 

  bagof(B-C, A^(ground_list(N,B),ground_list(N,C),foo(A,B,C)), BCS).



as(N, AS):- 

  bagof(A, B^C^(ground_list(N,B),ground_list(N,C),foo(A,B,C)), AS).



proof(N):-

  as(N,AS), bcs(N,BCS), 

  length(AS,N1), length(BCS, N2), N1 =:= N2.

This compares the number of successful B-C combinations overall with the number of As they produce. Equality means one-to-one correspondence.

And so we have,

2 ?- proof(2).

true.



3 ?- proof(3).

true.



4 ?- proof(4).

true.



5 ?- proof(5).

true.

And so for any N it holds. Getting slower and slower. A general, unlimited query is trivial to write, but the slowdown seems exponential.